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Structural features in a coalfield play an important role in the location of coal 
seams. In India, with the recent emphasis on power generation through thermal 
power plants, coalfields have been the centres of intense activity. Singrauli coalfield 
is taken as a case study as it has the unique distinction of having the thicJtest coal 
seam in India. This region has five super thermal power stations which act as nuclei 
for industrial development. The present work deals with digital image processing 
techniques for extracting the structural features in the coalfield. In addition to 
several imag(‘ processing techniques involving image smoothing, filtering and edge 
linking, a ndatively new concept of image transformation based on mathematical 
morphology has also been detaik'di A comparison of the performance of these 
techniques has been attempted and a final scheme of image processing for lineament 
detection is propo.sed. Results obtained have indicated a distinct improvement 
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over the available published work on the structural details within this coalfield. 
The methodology and the sequence of operations proposed are applicable to other 
coalfields as well. The softwares ased at different stages of work have been developed 
by the author during the course of this study. 

The work contained in the thesis has been organized into seven chapters. 

In Chapter 1, the topic has been introduced and an overview of the geology of 
the Singrauii coalfield is presented. The scope of the present work has also been 
indicated. In the north-eastern part of the coalfield taken up for study, the bound- 
ary between Lower Gondwana formations and the Precambrian metamorphic rocks 
is faulted. The Gondwana sequence is composed of the basal Talchir Formations 
(Lower Permian) followed by the Barakars, Barren Measures and Raniganj Forma- 
tions. The formations are essentially sandstones and shales. Coal seams in the area 
are as.sociated with the Baraloir and the Raniganj Formations. Several faults and 
lineaments are present in all these formations. Geological details along with the 
field features of all thes(> formations are presented in tliis chapter. 

Description of data nsed and details of the preprocessing techniques adopted 
are presented in Chapter 2. -Remotely sensed data firom two different satellites (IRS 
IB-LISS II and LANDSAT-TM) have been used. Preprocessing of satellite data 
has been carried out to establish compatibility as the spatial resolution, sensing 
geometry and system configuration differ. To reduce the problem associated with 
dimensionality of data, the information distributed in number of bands has been 
compressed into fewer number of components using Karhunen-Loeve transforma- 
tion. A base map of the study area is prepared using the satellite data as well as 
the available geological maps. 

Using the preprocessed IRS data, two false colour composites (FCC) have been 
prepared with combinations of bands 1 to 4. On the basis of the.se composites, the 
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lithological demarcation is carried out. Geomorphic features of the area are mapped 
on the basis of a FCC prepared using the principal components computed earlier 
(Chapter 2). A land use / land cover map of the study area has also been prepared 
which has brought out clearly the environmental impact of open cast coal mining. 
The details of all these aspects are presented in Chapter 3. 

Chapter 4 deals with digital enhancement of structural features. Image smooth- 
ing techniques have been applied to reduce noise in the satellite imagery. It has 
been realized that pseudomedian filtering offered desirable results. Edge enhance- 
ment has been carried out using derivative filtering of the first and second order, in 
addition to the entropy based filtering. Under the first order derivative filtering, or- 
thogonal and directional filtering are employed while Laplacian filtering is adopted 
under the second order filtering. Of all the operators used for edge detection, Prei- 
Chen operator of the orthogonal derivative filtering has yielded ber;t results. This 
image is thresholded and the output image is proces.sed using edge linking operation 
to render the final edge image. 

In Chapter 5, an edge detection approach based on mathematical morphology is 
discuSvSed. Image smoothing using neighborhood criterion has been performed for 
noise suppression. Transformations developed are dedicated to edge detection and 
involve the fundamental morphological operatioas such as dilation, erosion, opening 
and closing). Top hat transformations are developed based on opening, closing and 
image {liff<>rencing to enhance the edges in various directions. In the present work, 
it is seen that the top hat transformation used along with image superimposition 
technique provides a better edge image for the interpretation of lineaments. In the 
above transformations, several structuring elements of different shapes and sizes 
have been tested and the one yielding best result has been chosen. 

Chapter 6 deals with the lineament mapping and analysis. Based on the image 
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processing tools employed (Chapter 4 and 5) two lineament detection schemes - 
designated in the present work as Digitedge and Morphedge Schemes are utilized 
in preparing the lineament map. A statistical study of the characteristics of the 
lineaments has also been undertaken. The statistical analysis of lineaments involved 
tests for randomness, trend and specific orientation. Based on the qualitative and 
statistical analysis, the Morphedge Scheme is considered more appropriate for 
lineament mapping compared to the Digitedge Scheme. For delineation of faults, a 
scheme based on the morphological transformations has yielded encouraging results. 
Structural map of the area has been prepared using output obtained from this 
scheme. Several faults and lineaments not reported earlier could be delineated 
which were confirmed through limited field checks. 

Chapter 7 contains an overall siunmary of the present work. Based on the results 
obtained, suitability of various image processing techniques for structural mapping 
has been discussed. 
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Chapter 1 


INTRODUCTION 


1.1 General 

Among the three means of power generation, namely : the nuclear, hydro and 
thermal, in India the emphasis has largely been placed, of late, on the thermal power 
generation. In this context, Singrauli coalfield has gained pronodnence and presently 
occupies unique place in the coal map of India. It has the unique distinction of 
having the thickest coal seam known as the Jingurdah scam with a thickness of 
132 m. Coal is mined through open cast mining and the production is substantial 
(around 30 MT). This basin has also the unique distinction of having five Super 
Thermal Power Station.s and a host of industrial establishments being located with 
these power plants as the nuclei. Thns, the entire scenario in regard to the economic 
development of this region in India has gone through a sea change. In the last thrc^ 
decades geological mapping has been in progress and extensive drilling has enabled 
the delineation of several coal seams including the Jingurdah coal seam. For detailed 
exploration, this coal field was divided into eleven mining blocks and over the years, 
assistance from international agencies has also been availed for developing the coal 
mines. While conventional geological mapping and drilling have been the tools for 
the exploration, enough scope exists to examine the region using satellite data. Till 
date, except for scattered reports on visual interpretation of imagery, no work has 
been reported dealing with the digital analysis of satellite data for this coalfield. The 
study of remote sensing data facilitates an understanding of the regional picture. 
With the application of proper digital image processing techniques, the geological 
and structural details can be selectively enhanced. With this idea in view the 
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present work was undertaken. Satellite data as obtained through IR5-1B (LISS-II) 
and LANDSAT-TM pertaining to the year 1992 has been processed in the present 
work. The relevant software involved in the image processing and enhancement at 
different stages has been developed by the author. 

1.2 Singrauli Coal Basin - An Overview 

The Gondwana coal basins in Peninsular India occur as isolated intracontinental 
rift basins in the Precambrian terrain. These basins are irregular, elongate and oval 
with differing dimensions. They bear several similarities. Cyclic disposition in 
these basins is very characteristic. During deposition, basin subsidence has been 
the general feature. The basins are marked by prominent liigh angle normal faults 
that limit their boundaries (Basu and Srivastava, 1980). In Peninsular India, these 
basins occur along two preferred orientations namely EW and SE-NW (Fig. 1.1). 
Singrauli coal basin is an important member of this group. 

The basin has an area of 2200 sq. km of which approximately 4 percent is located 
in Sonbhadra district of U.P. with the rest extending into Sidhi and Sahdol districts 
of M.P. The basin is divided into two sub-basins, known as the Moher sub-basin 
in the East and the Main sub-basin in the West, the dividing line approximately 
conforming to the longitude 82° 30' E. In the Moher sub-basin, mining is in progress 
in the NE part. The coal field stands out as a prominent plateau in contrast to the 
adjacent plains. The coal field is well connected through transportation network 
with several places as Satna (307 km) in M.P. as well as with Pipri (40 km) , 
Renukoot (50 km), Mirzapur (193 km) and Varanasi (210 km) in U. P. It has railway 
links with Mirzapur and Katni. The climate is of tropical monsoon type with severe 
winter and extreme summer temperatures. Annual rainfall varies between 130 and 
150 cm. 
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Fig. 1. 1 -.Location of major Coalfields of India (Source: CMPDl, 1984) 



Geology of the study area 

As indicated in Fig. 1.2, the NE part of the basin is chosen as the study area. 
The sequence within the coal field in this region is essentially of Lower Gondwana 
(Upper Carboniferous - Permian) formatioas. An EW trending fault located in the 
northern part of the area limits the extension of these formations and marks the 
boundary with the Precambrian metamorphic rocks. The Precambrian formations 
are composed of gneisses, slates and phyllites. 

Within the gneissic terrain, patches of granite and pegmatite (Fig. 1.3) are 
also evidenced. In the area adjacent to the Bijul river (north of Singrauli Railway 
station), disturbance in the foliation has been evidenced (Fig. 1.4). In the vicinity of 
the boundary, the metamorphics outcrop in a series of parallel and elongate ridges. 
The general trend of foliation is along N 85 ° with high dip angles towards South 
(Fig. 1.5). Within the slates and phyllites, faulting has resulted in the dislocation 
of foliation trends and crushing of the rocks (Fig. 1.6). 

The generalized sequence (after Joshi and Pant, 1971) of the lower Gondwana 
formations in the study area is indicated in Table 1.1. Continued subsidence of the 
sediments along the boundary fault appears to have increased the tilt which has 
resulted in these formations dipping towards the fault. Accumulation of greater 
thickness of sediments close to the boundary faults in the Gondwana basins has 
been attributed to subsidence (Basu and Srivastava, 1980). 


4 



LEGEND 



Fig. 1.2: Geological mop of N-E Singrouli Coalfield (Source: Sinha,et al.,1982 ) 



Fig. 1.3 Pegmatite in gneissic terrain ( 5 km north of Morwa). 



Fig. 1.4 Brecciated gneiss in the fault zone. 
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Fig. 1.6 Disturbed trends in foliations of phyllites in the faulted zones. 
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Table 1.1 : Geological sequence of the Singrauli basin (after Joshi and 

Pant, 1971) 


Age 

Formation 

Lithology 

Thickness 

(m) 


Raniganj 

Sandstone, Shale 

380 



and Coal Seams 



Barren Measmes 

Coarse Sandstone 




Purewa Top Seam 

8-12 

Permian 


Sandstone 

0-60 


BaraJkar 

Purewa Bottom Seam 

10-14 



Sandstone and Shale 

45-75 



Kota Seams 

1-3 



Sandstone and Shale 

150-250 

Upper 

Talchir 

Shale and Sandstone 

not proved 

Carboniferous 





Talchir formations 

These are exposed mostly in the eastern and the south-eastern part of the study 
area. The Talchir sequence has a basal diamictite horizon composed of striated 
pebbles of sandstone, jasper, gneisses and basic rocks (Murthy,1957). The pebbles 
are of varying sizes and, at places boulders with striations indicating the glacial 
origin have also been reported (Pant and Shome, 1963). A thin horizon of diamictite 
evidenced in the southern part of the study area has elongated and rounded pebbles 
of differing sizes (Fig. 1.7). The Talchir sandstone is cross-bedded with current 
directions to NW (Raja Rao, 1983). The slate horizon in Talchir sequence is very 
thin. 
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Fig. 1.7(a) Pebble horizon in Talchir sandstone near Dheki. 
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Fig. 1.7(b) Pebbles in Talchir diamictite hoiiZOi near Jaitpur. 
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Fig. 1.8 Loosely cemented Barakar sandstones near Jayant block. 



Fig. 1.9 Barakax sandstones with ferruginous cement ( 2 km north of Singrauli) . 
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Fig. 1.10 Sandstones exhibiting slumping of beds and cross - bedding near 

Dudhichua. 



Fig. 1.11 Gray arenaceous shales overlying Barakar sandstones near Khadia. 


12 





Fig. 1.12 Coal seam underlying Barakar sandstones in Dudfaichua mines. 



Fig. 1.13 Formation of jhama due to burning of coal (Dudhichua). 
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Raniganj formations 

In the Jingurdah area, these, formations have an arcuate disposition and attain 
a maximum thickness of 400 m as reported by earlier workers. The sequence essen- 
tially is of sandstones associated with clays and carbonaceous shale. The sandstone, 
in general, is fine grained and has a massive appearance (Fig. 1.14). 

In the Jingurdah area, two coal seams occur within the Raniganj formations. 
The Lower seam, known as the Jingurdah bottom seam is reported to be of thick- 
ness upto 15 m while the upper seam designated as the Jingurdah top seam has a 
thickness of the order of 135 m. This is the single thickest coal seam reported so 
far in the coalfields of peninsular India. These two coal seams have a semi-elliptical 
garland shape and are offset by several cross cutting radial faults. 



Fig. 1.14 Fine grained comapct Raniganj sandstone with underlying coal seam 

(Jingurdah Mines). 
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For the development of the coal in the area, the Singrauii coal field has been 
divided into eleven mining blocks - designated as the Jingurdah, Gorbi, Kakri, 
Bina, Marrak, Khadia, Dudhichua, Jayant, Nigahi,Amlohri and Moher blocks. The 
reserves of coal in this coalfield axe reported to be aroimd 10851 MT (CMPDI, 
1984). The present study area encompasses nine out of these eleven blocks with the 
exclusion of Gorbi and Moher blocks. 

1.3 Scope and Objectives of the Present Study 

As already indicated, no published reports are available dealing with digital 
processing of satellite data for this coal basin. The study area offers excellent scope 
for investigations in this direction. The area exhibits distinct geomorphological 
features. The coal field area occurs as a plateau with scarp faces extending to the 
south of the area reflecting different stages of peneplanation. Erosional features 
like mesa and butte are evidenced at several places in Barakar formations. The 
Precambrian metamorphic formations stand out as prominent linear ridges in the 
northern part of the area. In addition to the EW boundary fault, several faults are 
present in the Gondwana formations. Some of them, as in the Jingurdah coal mines 
are radial in nature and offset the coal seams. The lineaments characteristic within 
the coal field area are different compared to the adjoining metamorphic terrain. 
Selective enhancement techniques are thus called for to delineate all these features 
and to prepare a detailed structural map. With this in view, the following aspects 
are taken up for detailed study leading to structural mapping of coal field area: 

1. Evaluation of different image smoothing, edge detection and edge lining tech- 
niques for enhancement of structural features and to finally evolve a scheme 
of operation. 
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2. Application of mathematical morphology based image processing techniques 
involving fimdamental operations of dilation, erosion, opening and closing for 
delineation of lineaments. 

3. To evaluate the above two image processing schemes in terms of their better 
detection of structural features and computational efficiency. 

4. Preparation of a structural map of the study area using the edge image re- 
ceived from these schemes. 

In the present work, structural mapping is confined to lineaments including 
faults. 


1.4 Organization of the Work 

The material contained in the thesis is organized and presented under seven 
chapters. Chapter 1 deals with the introduction. Details of the Singrauli coal basin 
and study area are given. The scope and objectives of the study are defined. In 
Chapter 2, the data used and the preprocessing techniques adopted are presented. 

Based on the false coloured composites created with different band combinations 
and principal components, geological and geomorphological as well as land use maps 
are prepared. The details of the same are presented in Chapter 3. Chapter 4 
deals with digital enhancement of structural features using image smoothing, edge 
detection and edge linking techniques. 

In Chapter 5, an edge detection approach, based on mathematical morphology 
is presented. An evaluation of different operations involved, is attempted. Based 
on the techniques developed in the previous chapters, two schemes of image pro- 
cessing for lineament mapping are proposed and output images have been used for ; 
preparation of maps. A statistical analysis of lineament has been carried out for 
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the area. The details are presented in Chapter 6. The findings of the present work 
are s umm arized in Chapter 7. 


17 



Chapter 2 


PREPROCESSING OF DATA 


2.1 Introduction 

During the last two decades, launching of several unmanned satellites into the 
orbit around the Earth has resulted in obtaining the satellite data dealing with fea- 
tures on Earth’s surface. Information is now available from several sensor platforms 
such as LANDSAT - TM, SPOT, ERS and IRS - IB. In the present work, multi- 
band and multiresolution data have been selected from two sensor systems namely 
the LANDSAT - TM and IRS - IB. These data vary in their spectral ranges and 
spatial resolution. The orbital path and elevation of the platform is also different in 
each case. In addition, the data have been recorded on two different dates. Hence 
to integrate this dataset for further processing, geometric rectification and image 
resampling were performed. 

To take care of the dimensionality of multichaimel dataset, a transformation 
based on principal component analysis has also been performed. This trans- 
formation has reduced the number of spectral bands and provided the maximum 
details of the Earth cover types. All these preprocessing techniques adopted in the 
present work are detailed in this chapter. 

2.2 Data Acquisition 
Satellite data 

In the present work, data from LANDSAT - 5 and IRS - IB were used. The 
LANDSAT-5 launched on March 1, 1984 and moving in a Sun synchronous orbit 
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at an altitude of 705 kms with 16 days repetitive coverage having a field view of 
14.93° carries a four - channel naultispectral scanner and a thematic mapper (TM). 
The TM is a scanning optical - mechanical sensor system that records reflected 
and emitted energy in the visible, reflective-infrared, middle-infrared and thermal- 
infrared regions of the electromagnetic spectrum and the spatial resolution is 30 
m for Bands 1 to 5 and 7 while for Band 6 it is 120 m. The data collected has 
higher spatial, spectral and radiometric resolution than the one from MSS. Details 
regarding LANDSAT - TM spectral bands and their applications are presented in 
Table 2.1. 

On the coverage index map, the study area is located by a path number - row 
number 142 - 043 and the data pertains to April 21, 1992. 

Indian Remote Sensing Satellite (IRS - IB) was launched on August 29, 1991 
in a Sun synchronous orbit at an altitude of 904 kms. This satellite carried one 
pushbroom camera (LISS - I) of 72.5 m resolution and two cameras (LISS - II) of 
36.25 m resolution and images obtained are of four spectral bands in the visible and 
near IR region. The spectral characteristics of IRS - IB and SPOT are presented 
in Tables 2.2 and 2.3. 

The study area is covered on IRS - IB index map by path number-row nximber 
23 - 51 and the data pertains to April 25, 1992. The data from IRS - IB covers 
essentially the same spectral range as that from SPOT. Although the resolution of 
SPOT - MLA data is better compared to IRS - IB LISS - II, the data from SPOT 
is either of one band (SPOT - PLA) or of three bands (SPOT - MLA) while that 
from IRS - IB is of four bands covering the same spectral range. Thus IRS - IB 
has an advantage in providing more details. 
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Table 2.1 LANDSAT - TM spectral bands and principal applications 
(Jensen, 1986) 


BANDS 

SPECTRAL 

RANGE 

(micron) 

SPECTRAL 

LOCATION 

PRINCIPAL APPLICATION 

1 

0.45 - 0.52 

BLUE 

Sensitivity to sedimentation studies, 

coastal water mapping. 

soil / vegetation discrimination, 

forest type mapping. 

2 

0.52 - 0.60 

GREEN 

Green reflectance of healthy vegetation, 

vigor assessment in vegetation. 

3 

0.63 - 0.69 

RED 

Sensitivity to chlorophyll absorption by 

vegetation, differentiation of soil 

and geological boundaries. 

4 

0.76 - 0.90 

NEAR 

INFRARED 

Sensitivity to green biomass 

and moisture in vegetation delineating the 

water bodies and soil moisture mapping. 

5 

1.55 - 1.75 

MID - 

INFRARED 

Useful for vegetation moisture 

content and soil moisture , 

discrimination of snow from cloud. ' 

6 

10.4 - 12.5 

THERMAL 

INFRARED 

Useful in vegetation moisture analysis, 

soil moisture discrimination, 

and thermal mapping applications. 

7 

2.08 - 2.35 

MID - 

INFRARED 

Useful for discrimination of minerals 

and rock types. Also sensitive 

to vegetation moisture content. 
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Table 2.2 IRS spectral bands and principal applications (Dept, of 
Space, India, 1989) 


BANDS 

SPECTRAL RANGE 

(microns) 

SPECTRAL 

LOCATION 

PRINCIPAL APPLICATION 

1 

0.45 - 0.52 

BLUE 

Sensitivity to sedimentation, 

deciduous / coniferous forest 

cover discrimination. 

2 

0.52 - 0.59 

GREEN 

Green reflectance of healthy 

vegetation. 

3 

0.62 - 0.68 

RED 

Sensitivity to chlorophyll 

absorption by vegetation, 

differentiation of soil 

and geological boundaries. 

4 

0.77 - 0.86 

NEAR- 

INFRARED 

Sensitivity to green biomass 

and moisture in vegetation 
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Table 2.3 SPOT - MLA spectral bands and principal applications 


(Barret and Curtis, 1992) 


BANDS 

SPECTRAL RANGE 

(|im) 

SPECTRAL 

LOCATION 

PRINCIPAL APPLICATION 

1 

0.50-0.59 

GREEN 

Sensitivity to sedimentation, 

deciduous / coniferous forest cover 

discrimination. 

2 

0.61-0.68 

RED 

Reflectance of healthy vegetation. 

3 

0.7()-().89 

RliFLBCTED 

-INFRARED 

Sensitivity to chlorophyll absorption 

by vegetation, dil ferenlialion of soil 

and geological boundaries. 



Table 2.2 IRS spectral bands and principal applications (Dept, of 
Space, India, 1989) 


BANDS 

SPECTRAL RANGE 

SPECTRAL 

PRINCIPAL APPLICATION 


(microns) 

LOCATION 


1 

0.45 - 0.52 

BLUE 

Sensitivity to sedimentation, 

deciduous / coniferous forest 

cover discrimination. 

2 

0.52 - 0.59 

GREEN 

Green reflectance of healthy 




vegetation. 

3 

0.62 - 0.68 

RED 

[ Sensitivity to chlorophyll 




absorption by vegetation, 




differentiation of soil 




and geological boundaries. 

4 

0.77 - 0.86 

NEAR- 

Sensitivity to green biomass 



INFRARED 

and moisture in vegetation 
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2.3 Data Retrieval from CCT 

The data from NRSA was supplied on the computer compatible tapes (CCT). 
It was read using a computer software on MICRO VAX System. Since the data on 
CCT is written in a particular standard format (Ramkrishnan, 1988), a computer 
software in VMS FORTRAN has been prepared by the author based on the specified 
format. There are five files on the CCT organized in a superstructure format as 
described below: 

• Volume Directory 

• Leader 

• Image Data 

• Trailer 

• Null Volume 

The volume directory file contains information relating to the tape data. It 
is composed of volume description record, a series of file pointer records and a 
text record. The leader file constitutes file description record and three types of 
data record. The data types are header, ancillary and annotation. The image file 
consisting of file description record and actual image data record is followed by a 
trailer record which finally ends with a NULL volume directory file. 

The data in the image file is written in two formats either band interleaved by 
line (BIL) or band sequential (BSQ). In present case the data was procured in BIL 
format. The first record on the image data file is the file description record of 2520 
bytes. In addition, there are 10,000 image data records each of 2520 bytes. The 
main data on each of these image records are contained in bytes 33 to 2382. A 
layout of user CCT is shown in the Fig. 2.1. 
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Figure 2.1 : Layout of user CCT of BIL format 
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A computer software to retrieve the data of different bands of a specified area or 
subscene as written by the author is presented in the Appendix - A. With the help 
of this software, data on CCT has been accessed and the exact portion covering the 
study area has been windowed and the data read and recorded on to the floppies 
for further processing on a personal computer. The computer system used is of PC 
386 with super VGA card and a 19" colour monitor. 

The satellite data of the study area comprising of four bands of IRS and one 
Band of LANDSAT - TM used for digital image processing and analysis are shown 
in the Figs. 2.2 to 2.6 . 

2.4 Preprocessing of Satellite Digital Data 

Data supplied to a user is rectified with respect to systematic and non - sys- 
tematic errors. These errors are rectified before the data being supplied to the 
user. However, the data used belongs to two satellites, the present problem is to 
eliminate the systematic and non - systematic incompatibility. The systematic in- 
compatibility is attributed to platform velocity, Earth rotation, scanning system, 
spatial resolution and spectral resolution. The non - systematic incompatibility is 
due to altitude, attitude and atmospheric factors. 

As already stated that satellite data being used two sensor systems is having 
different spatial resolution. The data from these sensors are of different dates. 
The platforms of these sensors are placed at different altitudes and orbits thus 
involving atmospheric variability. The data can be integrated only after performing 
corrections to reduce the incompatibility. For this purpose geometric rectification 
and image resampling have been performed. 
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Fig. 2.3 IRS IB LISS II image o idy area, Band l. 
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Fig. 2.5 IRS IB LISS II image of study area, Band 4 
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The LANDSAT - TM data (with a resolution of 30 m) is to be transformed into 
the grid size of IRS - IB data (with a resolution of 36.25 m). For this purpose the 
geometric transformation is carried out as follows : 


■^M = ^ [/ (^TM, -RiRSi CiRs)] 

(2.1) 

C'TM=^^[/(<^TM,RlRS,aRs)] 

(2.2) 


where 

Rirs is line number of IRS - IB image, 

CiRs is pixel number of IRS - IB image, 

Rxm is line number of LANDSAT - TM image, 

Ctm is pixel number of LANDSAT - TM image, 

is line number of geometrically corrected LANDSAT - TM image and 
Ctm i® pixel number of geometrically corrected LANDSAT - TM image. 

For obtaining a transformation of this t3q5e, coordinates of ground control points 
(GCPs) present on both the imagery are required. Eight ground control points were 
identified both on LANDSAT - TM image and on IRS - IB LISS - 11 image. For 
the geometric correction, image of any band of IRS can be taken. In the present 
operation Band 3 image of IRS - IB (Fig. 2.4) was utilized as the same is used 
subsequently. LANDSAT - TM image utilized for transformation to IRS - IB 
grid domain is presented in Fig. 2.6. The size of Band 3 IRS image is 512 x 512 
whereas the size of uncorrected LANDSAT - TM image is 632 x 616. Three possible 
approaches are available for establishing a mapping relationship between line and 
pixel numbers of LANDSAT - TM and IRS - IB imagery. They are generally linear, 
quadratic or cubic. The ground control points required for these three cases should 
be greater than 3, 6 and 10 respectively. For greater accuracy either quadratic 
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or cubic relationship is preferred over the linear type. In the present work, the 
quadratic relation was adopted since the available ground control points were limited 
within the study area. The quadratic mapping function is of the following form : 

Rtm — K + b\CiRs + h^Rms + hCjjis + biCjjtsRiRS + hRjjts ( 2 - 3 ) 

and, 

Ctm = ^>0 + b'lCiRs + b2RiRs + K^irs + b'^CjRsRiRs + KR]rs (2-4) 

Eight GCPs from the LANDSAT - TM and IRS - IB images were picked up 
(Table 2.4) on the basis of which the coefficients of the quadratic mapping func- 
tion (Table 2.5) for transforming LANDSAT - TM image coordinates to IRS - IB 
coordinate system were computed. 

After the geometric rectification of the LANDSAT -- TM image into IRS - IB 
grid system, the resampling of the LANDSAT - TM image is performed. 

Table 2.4 : Ground control points selected for the study area 


Point 

number 

Line number 

(LANDSAT) 

Pixel number 

(LANDSAT) 

Line number 

(IRS) 

Pixel number 

(IRS) 

1 

573 

449 

433 

332 

2 

139 

486 

70 

357 

3 

298 

579 

205 

439 

4 

406 

578 

296 

439 

5 

184 

325 

105 

220 

6 

189 

355 

112 

245 

7 

386 

464 

276 

343 

8 

561 

247 

420 

161 
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Table 2.5 : Coefficient of quadratic mapping function for geometric correction 


Coefficient 

Ctm 

coefficient 

Rtm 

bo 

53.7774 

b'o 

80.2012 

i.. 

1.2139 

b\ 

-0.0658 

i>2 

0.0234 

b'2 

1.1259 

h 

-0.0001 

b'z 

0.0001 

64 

0.0 

b\ 

0.0001 

65 

-0.0001 

b'. 

0.0001 


2.4.2 Image resampling 

Once the transformation of grid system has been achieved, it is necessary to 
read the corresponding grey level values and put them in the new grid domain. 
The original grid as well as the geometrically corrected grid are shown in Fig. 
2.7. As indicated in the figure, the original grey level values can not be directly 
associated to any of the corrected grid cell. For properly assigning the appropriate 
grey level values in the new grid, the following techniques are available (Lillesand 
and Kiefer, 1979). 

1. Nearest neighbourhood technique 

2. Bilinear interpolation technique 

3. Cubic convolution technique 

The nearest neighbourhood technique has been applied in present case. In 
this method, the grey level of the nearest cell is adopted for the corrected grid cell. 
Therefore in the example shown in Fig. 2.7, the grey level of cell 3-3 will be 
recorded in the shaded cell of corrected grid. Here the grey levels are not altered. 
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Fig. 2 . 7 : Nearest neighborhood process for resampling. 




T his technique is computationally very simple and straightforward. The bilinear 
interpolation and cubic convolution techniques require large computational time 
compared to this operation. The grey levels are also modified. Because of t his , 
image may look little smoother but with the transformation of grey level the basic 
aim behind selecting a particular spectral band is denied. The LANDSAT - TM 
Band 7 is selected for study as it is very useful for geological feature interpretation. 
If there is any kind of modification in these spectral responses, then the usefulness 
of the image for geological applications may severely be affected. Thus, the process 
of nearest neighbourhood sampling is well justified in present application. 

Using this method, the LANDSAT - TM image is transformed in to identical 
grid as of IRS - IB LISS - II and the output LANDSAT - TM image is shown 
in Fig. 2.8. This LANDSAT - TM image has been made use of in subsequent 
interpretation and processing along with the four IRS -- IB LISS -- II bands. 

2.5 Data Compression 

In the present case the main features targeted for identification are geological. It 
is also understood that vegetation, drainage patterns and moisture content of surface 
contribute significantly in the identification of geological details. Incorporating all 
these details in the best form is not possible by selecting a single spectral band. 
At the same time it is very diflScult to handle large data simultaneously. Another 
problem of concern, is that the image objects often have intermixing nature of grey 
level values. This problem is severe in a highly correlated dataset and imposes a 
problem of misclassification and misinterpretation. To maximize the variance and 
minimize the correlation, the technique of principal component analysis (Jolliffe, 
1986; Jensen, 1986) was utilized. The principal components are computed using 
the Karhunen - Loeve transformation (Moik, 1987). 
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The main objective of principal component analysis is to transform a highly 
correlated n - dimensional data set into a lesser dimension feature space by redis- 
tributing the variance and information content. Number of bands or dimension of 
data should be equal or more than the desired components. 

For understanding the concept of principal component analysis a two - dimen- 
sional case has been explained that can also be extended to the multi - dimensional 
case. For two bands bi and the principal components pi and p 2 can be computed 
as 

pi = Cubi ■+C 12&2 (2-5) 

P 2 = C 2161 + 022^2 (2-6) 

This can be expressed in the matrix form a 5 , 

Pi _ cii C12 bi 

_ P 2 J [ C 21 C 22 J [52 

By a graphical representation this concept can be better appreciated as shown 
in Fig. 2.9. The data of two spectral bands of an image are plotted on two axes bi 
and 62 (Fig. 2.9 a). The scatter plot is reflecting a highly correlated dataset and 
hence the image is of very little information value. To augment the information 
value of an image, the correlation of dataset has been reduced by transformation to 
maxinoize the variance. In the two - dimensional case, first a linear transformation 
is applied by linear translation of the origin of previous coordinate system. For 
this purpose mean values of bands, bi and 62 are determined as pi and p 2 - A new 
coordinate system based on b\ and 62 is then derived by a simple expression 

5 j —bi+pi 

62 = 52+ P2 
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(a) Scatter plot of Band 1 and 2; (b) Linear translation 

of the origin and axes; (c) Rotation of axes to map 
the data in the direction of maximum variance. 



This is illustrated in Fig. 2.9 (b). In order to obtain the maximum variance 
along one principal component, the newly found axis is again rotated by an angle 
0 so as to align b[ along the maximum variance direction which is along the major 
axes direction in the present figure. This component is the first principal component 
(PCi) and the component orthogonal to it is the second principal component (PC2). 
The rotation transformation is illustrated in Fig. 2.9 (c). In order to reproject 
the original dataset onto the new {PCi) and {PC2) axes, certain transformation 
coefficients are needed which are derived from the covariance matrix of the original 
dataset. 

Similarly on a multichannel dataset which is highly correlated, a transformation 
has to be applied to produce a new multichannel uncorrelated dataset having an 
ordered set of variance. The transformation for a multichannel spectral case is 
carried out in the following steps : 

(1) The mean and standard deviation values for the given spectral bands TM 
(Band 7 ), IRS (Band 1), IRS (Band 2), IRS (Band 3 ) and IRS (Band 4 ) are deter- 
mined. 

( 2 ) Variance - Covariance matrix for this dataset is computed (Table 2 . 6 ). 

( 3 ) The correlation matrix is computed (Table 2 . 6 ). 

( 4 ) The eigen values and eigen vectors are calculated such that 


[EV][COV][EV'^] = 


An 0 0 0 0 

0 A22 0 0 0 

0 0 A33 0 0 


0 0 0 A44 0 

0 0 0 0 A55 

where, A„„, n = 1 . . . 5 are eigen values for Bands 1 to 5 . Here, Band 1 designated 


( 2 . 8 ) 
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for TM Band 7, and Bands 2 to 5 are for Bands 1 to 4 of IRS. The right hand side 
matrix is the diagonal matrix having respective eigen values as diagonal elements. 
EV is eigen vector. 

EV = K,] 

where, r = 1 to 5 Bands and s = 1 to 5 components. 

EV^ is the transpose of matrix EV 

Percent of total variance (%pc) contained in each principal component is deter- 
mined by the expression 

%pc = 100 2.9 

where Apc is the eigen value for concerned principal component. In present case n 
= 5. 

By evaluating the correlation of each band with each component it is possible to 
determine how each band is associated with each principal component or its factor 
loading. 


fr. = ( 2 . 10 ) 

<^r 

The computed percent variance, eigen values, eigen vectors and factor loadings 
are listed in the Table 2.7. 

In above case five principal components have been computed for the five spectral 
channels in consideration. As evident firom the Table 2.7, the maximum information 
(upto 61.64 percent) is contained in the first principal component. The image of 
this principal component is shown in the Fig. 2.10. The first and second compo- 
nent together offer information upto 92.80 percent. In addition, the fourth principal 
component is integrated with the first two. Although the fourth principal compo- 
nent exhibits only 1.36 percent variance, its image (Fig. 2.11) depicts a detailed 
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lineament fabric of the study area. The first, second and fourth principal compo- 
nents were used in this study to generate FCC which has been used in the study 
for further interpretation. 

Table 2.6 ; Variance - Covariance matrix and Correlation matrix of 
original dataset 


Variance — Covariance Matrix 



TM7 

IRSl 

IRS2 

IRS3 

IRS4 

TM7 

481.35 

-31.17 

-20.52 

16.22 

76.26 

IRSl 

-31.17 

56.38 

48.73 

70.31 

20.15 

IRS2 

-20.52 

48.73 

48.09 

75.06 

24.83 

IRS3 

16.22 

70.31 

75.06 

138.31 

65.41 

IRS4 

76.26 

20.15 

24.83 

65.41 

85.70 


Correlation Matrix 



TM7 

IRSl 

IRS2 

IRS3 

IRS4 

TM7 

1.0000 

-0.1892 

-0.1349 

0.0629 

0.3755 

IRSl 

-0.1892 

1.0000 

0.9359 

0.7962 

0.2899 

IRS2 

-0.1349 

0.9359 

1.0000 

0.9204 

0.3868 

IRS3 

0.0629 

0.7962 

0.9204 

1.0000 

0.6008 

IRS4 

0.3755 

0.2899 

0.3868 

0.6008 

1.0000 
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Chapter 3 


GEOLOGICAL AND LAND USE MAPPING USING 

SATELLITE DATA 


3.1 Introduction 

For the demarcation of lithological boundaries, geomorphological features and 
structural details, false colour composites find extensive use. IRS data earlier has 
been used for similar purpose, effectively. Krishnamurthy (1990) used false colour 
composites generated from Bands 1, 2, 3 and 4 of IRS for discriminating litho units 
and for the identification of major geological and geomorphic features. Ganesh Raj 
and Gupta (1990) have prepared a hydrogeomorphological map based on FCC gen- 
erated firom IRS-IA LISS-II data. Newton and Boyle (1993) have utilized Landsat 
TM data for preparation of FCC of principal components to retrieve geological in- 
formation. According to them, such FCC is more useful for differentiating rock 
t 3 rpes compared to the FCC of SPOT data. 

To have an overview of the study area, false colour composites have been gener- 
ated and used in this chapter. These have been prepared using various combinations 
of foinr bands of IRS data. In addition, principal components calculated earlier us- 
ing four bands of IRS and Band 7 of LANDSAT TM, have been integrated to create 
FCC which has been utilized for geomorphological and land use mapping. 

3.2 Lithological and Structural Mapping 

False colour composites using IRS Bands 1, 3 and 4 (Fig. 3.1) and Bands 1, 2 
and 4 (Fig. 3.2) ha^been generated on VAX 11/780 computer system. 

A 
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Fig. 3.1 False colotir composite of IRS (Bands 1, 3 and 4) 
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Fig. 3.2 False colour composite of IRS (Bands 1, 


2 and 4) 
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These have been used along with Band 7 image of Landsat TM for the visual 
interpretation and mapping of geological and structural details (Fig. 3.3). In the 
study area, the different sedimentary formations of the lower Gondwana group 
could be clearly recognized and the boundaries marked. In addition, within the 
Precambrian metamorphic terrain, gneisses could easily be distinguished from the 
slates and phy Hites. 

The major EW trending boundary fault is seen clearly on the FCC along with 
another fault in metamorphic terrain in the eastern part of the area. Several faults in 
the Jingurdah region and in the SW part of the study area could be identified. Fault 
delineation in other parts of the area was not possible on FCC due to tonal contrast 
limitations. Methodology, presented in the subsequent chapters, has enabled a clear 
delineation of all the linear features in the entire study area. 

3.3 Geomorphological Mapping 

Classification of the area in terms of various geomorphic units is carried out. For 
this purpose, FCC (Fig. 3.4) is generated as already indicated in Section 3.1. The 
structural ridges in the region of slates and phyllites (for example: A2, B2,B1,C1) 
and in the gneissic terrain (D4,G4,H4,C4) are clearly visible. Drainage network of 
Bijul river in the NW part of the area (A4 to FI) with its dendritic pattern has 
been observed to flow at several places along the lineaments indicating a distinct 
structural control of drainage in the region (B3-C3, A4-B4, C3, D2-E2, El and Fl). 

As already indicated, the coal field area stands out, in contrast, as a plateau 
compared to the adjoining metamorphic region. The drainage pattern of Bijul river 
(B6-B4 and E5-E3) is in conformity with the same. Another drainage network Ballia 
nala (ElO, FlO, Ell, Hll, 111) flowing through coal mine areas of Dudhichua (E9) 
could be easily traced upto its final merger into the G. B. Pant reservoir. 
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The flow direction confirms the elevation gradient from Baraiar terrain to Talchir 
terrain. 

Pediment zones in several places (A5, B7, C6, E5, and E7) and pediplain areas 
(D6-E8) within Barakar terrain axe identified. Several Mesa and Butte features 
(A7, A9, B8, H5, H6, H7) are also mapped. The geomorphological map prepared 
is presented in Fig. 3.5. During field checks, typical geomorphic features observed 
in the field are presented in Figs. 3.6 to 3.10. 

3.4 Land Use Mapping 

The FCC (Fig. 3.4) of principal components has also been visually interpreted 
to prepare land use map (Fig. 3.11). Mining activity in the region has a significant 
impact on the ecology. As the coal mining is through open cast method, serious 
environmental hazards have resulted. These include depletion of vegetation or forest 
cover, degradation of agricultural land, piling of huge mine dumps (Fig. 3.12) and 
acid mine dradnage (Fig. 3.13) that pollutes the G. B. Pant reservoir. 





Fig* 3-5 Geomorphological map of the study area prepared from 
folse colour composite generated using principal 
components (Fig. 34) 






Fig. 3.6 Structural ridges in slates and phyllites near Khirba. 



Fig. 3.7 Structural ridges in gneissic terrain near Bijul river (North of Singrauli) 
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Fig. 3.8 Typical pediplain in the Barakar terrain near Merhauli. 



Fig. 3.9 A typical mesa in Barakar terrain near Nigahi. 
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Fig. 3.12 A mine dump in the Jayant block. 



B. Pant reservoir. 
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Chapter 4 


DIGITAL ENHANCEMENT OF LINEAMENTS 


4.1 Introduction 

Structural features such as lineaments, faults, joints, fractures, folds are man- 
ifested on the surface of the Earth by their specific signatures. Satellite images 
offer a global view of the Earth’s surface and facilitate the identification of these 
geological structural features. Due to various Earth cover types like forest, agri- 
culture and urban sprawl, it is often not possible to delineate structural patterns 
efficiently. Image processing techniques play an important role in enhancing the 
structural features in such cases. 

Satellite imagery is utilized for digital enhancement using filtering in frequency 
or spatial domain. Structural features are enhanced by enhancing the edges. Edges 
are defined as the boundary between two regions with relatively distinct grey level 
properties (Gonzalez and Wintz, 1977). The objective of filtering in present case is 
to enhance the edges. Choice of individual filters depends on specific requirement. 
If a priori knowledge of the Uneament pattern in the area is available, a filtering 
scheme based on this criterion may be adopted for the best result. Folds have special 
curvilinear nature while faults, fractures and joints approximately follow a straight 
line pattern. 

Edge detection is a process of applying a transformation {T) to an original image 
f{x, y) to obtain the output image g(x, y) containing edge informations. Mathemat- 
ically, It may be expressed as 


g{x,y) =r[f{x,y)] 


(4.1) 
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and hides important details. These techniques are computationally costly. The 
local operators are applied to the noisy images in a small neighbourhood of the 
concerned pixel in the computation. They have great potential of real time or 
near-real time implementation, because several pixels can be processed in parallel 
without waiting for their neighbouring pixels to be processed. Local processing 
can be implemented with better control over the preservation of subtle details and 
edges. Due to their computational efficiency as well as better result, local operators 
are gaining popularity over global smoothing (Lee, 1983). 

The most popular among local smoothing methods is the moving average method 
which is mathematically simple and useful in many applications. Median filter 
in one or two dimensions has become popular due to its characteristics of edge 
preservation and removal of spikes (Huang, 1981). The other popular smoothing 
operations used by many researchers are edge preserving smoothing of Nagao and 
Matsuyama (1979), the gradient inverse weighting scheme of Wang et al. (1981), 
the box-filtering algorithm (McDonnel, 1981) and the local statistics methods of Lee 
(1980, 1981). It is not appropriate to generalize their ranking since an algorithm 
effective for a class of images, may be ineffective for others. In the present study 
four image smoothing techniques namely moving average filtering, sigma filtering, 
median filtering and pseudomedian filtering are being discussed and compared^based 
on computational requirement as well as suitability. 

4.2.1 Moving average filtering 

An image /(x, y) is transformed to a smoothed output image /(x, y) whose grey 
level at each pixel location (x, y) is found by averaging the grey level values of the 
pixels of image contained in a predefined neighbourhood of central pixel at (x, y). 
The neighbourhood size varies with the nature of the image. The transformation is 
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obtained as 


= E fi^^y) (4-2) 

(i,y)es 

for X, y = 0, 1, . . . , {N — 1). S' is a set of coordinates of points in the neighbour- 
hood of the point (x, y) including (x, y) itself and N x N is the total number of 
points in the neighbourhood (Gonzalez and Wintz, 1977). It has been reported by 
Gonzalez and Wintz (1977) and Lee (1980) that this filter is useful in some cases 
but fails in many other applications due to its main drawback of blurring the edges. 

For a band 3 image of IRS, this method has yielded a better smoothed output 
over a 3 X 3 size window compared to a 5 x 5 size window. The image received in 
latter case is oversmoothed and appears blurred. The degree of blurring is found to 
be strongly proportional to the size of the neighbourhood used. 

Moving average method due to its drawback of blurring the edges, is unlikely to 
yield desirable smoothing output. 

4.2.2 Sigma jSltering 

Sigma filtering has an improvement over the moving average filtering. As already 
discussed, in the moving average filters the edges get blurred as it replaces the 
central pixel of the window by the neighbourhood average. This is the result of 
blind averaging of pixels which are not exactly the real members of a feature family. 
Nagao and Matsuyama (1979) and Lee (1981) have produced solutions to rectify 
this problem. 

Sigma filtering, also on sinoilar lines, uses the basic concept of the sigma proba- 
bility of Gaussian distribution. The primary aim in this filtering is to replace the 
central pixel of the kernel by the average of pixels which are members of the same 
feature class. This is attained by locating those neighbouring pixels having their 
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intensity within a fixed signaa range of the central pixel of the kernel. Replacing the 
central pixel of the kernel by the average of selected neighbouring pixels has also 
been suggested by various other researchers. The filter proposed by Nagao and Mat- 
suyauna (1979) replaces the central pixel of the kernel by the average of a subregion 
which has minimum variance. Lee (1981) has evolved a concept of using gradient 
information in his modified local statistics of neighbouring pixels. In the extended 
box-filtering algorithm the average is computed only for those neighbouring pixels 
falling within a fixed grey level domain (Pratt, 1991). Box-filter and sigma filter 
differ only in one way in that the box-filter has the grey level range fixed throughout 
the whole image while the sigma filter has freedom of adaptable grey level range. 

Several advantages of sigma filter are 

1. Ekige preservation : The noise in the image space and near the boundaries 
and edges will be smoothed without blurring the edge because only pixels 
falling on one side of the edges are utihzed while computing the average. 
Since not all pixels are included in average the blurring effect is avoided. 

2. Preservation of shape ; Sigma filtering does not resort to any biases for 
direction as adopted by Nageio and Matsuyama (1979) and Lee (1981) in their 
methodology. This helps in retaining the original shape of inherent details 
which is very important as for as the remote sensing imagery is concerned as 
a change in morphology affects the process of feature identification. 

3. Preservation of salient details : By selecting a proper window size, the 
salient subtle details of several pixel group pertaining to similar class and 
linear features of desired dimension can be designed to be preserved since 
only relevant pixels are included in computation of the average and not the 
background information. 
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4. Efficiency of computation : It is computationally more efficient since only 
simple comparison and addition operations are performed. 

The class of pixels in a very tiny neighbourhood, uncorrelated with the sur- 
rounding image space, is generally considered a noise. Most of the image noises 
are characterized by Gaussian distribution (Lee, 1983). The two-sigma probability 
for a feature space is defined as the probability of a random variable being within 
two standard deviations of its class mean. Here one can recall that the two-sigma 
probability for a one-dimensional Gaussian distribution is 0.955. This means that 
95.5% of random pixels lie within the range of two standard deviations. Any pixel 
falling outside the two-sigma range most likely comes from a different population 
and hence it should be excluded from the computation of the average. It is assumed 
that a pri,ori mean comes fix)m the two-sigma range of the grey level values of the 
pixels to be smoothed. 

Assuming that the noise is additive in nature with zero mean and standarad 
deviation of a, the sigma filtering (Lee, 1983) on a discrete image is performed as 
follows ; 

1. A grey level range (r^ + A, Vij - A) is formed, where n, is the grey level value 
of the pixel (i, j) and A = 2a. 

2. Compare the neighbourhood pixels in the specified window. 

3. Compute the sum of all the pixel lying in the specified grey level range. 

4. Compute the average by dividing the sum by the number of pixels included 
in the computation of the sum. 

5. Now fij — average. 
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6. Replace the kernel by fy. 


Mathematically the algorithm may be explained as follow : 

1, if (ry - A) < rki < (ry- + A); 
0, otherwise. 




Then, 


■ , — ( 4 - 3 ) 

By assuming a two-sigma range it is ensured to include 95.5% of the pixels 
from same feature space. At the same time it is also ensured not to include pixels 
representing high contrast edges such as spikes. Structural information from satellite 
imagery which is linear in nature, is mostly retained if it is two or more pixels wide, 
because only those pixels falling within the specified two-sigma range are included 
in the average. This is the main drawback as the scattered noise of one or two pixels 
size will not be smoothed, resulting in an undesirable output specially in case of 
noisy image. This problem is resolved by replacing the two-sigma average with the 
window average. If the number of pixels within the specified range, is less than a 
threshold value m then the step (5) is replaced by 


{ two-sigma average, if n > mj 

(4.4) 

immediate neighbour average, if n < m 
where, n is number of pixels in the specified range. 

The value of m should be carefully chosen to remove isolated spot noise without 

affecting linear features. For a 7 x 7 window, m should be less than 4 while for 

a 5 X 5 window it should be less than 3. It is found that subtle textures within 

two sigma-range are wiped out after few iterations (Lee, 1983). If conservation of 
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textural information is required which is vital in the case of remote sensing, a small 
A range and one or two iterations should be used. 

Generally local smoothing operators do not require any prespecified threshold 
or any parameter. Sigma filter requires specification of the grey level range and the 
size of the window. 

The computational ease and efficiency in addition to edge preserving nature 
and other advantages over contemporary smoothing algorithms makes sigma filter 
a distinct choice to be adopted in the case of remotely sensed imagery. The main 
drawback found in the present application is of selecting a common threshold m. 
The reason being that remote sensing image contains a wealth of information i.e. 
local variance in a neighbourhood changes abruptly. This results in turning the 
sigma filter into moving average filter most of the time. Thus all the advantages 
discussed above are not realized in the case of an overall image. 

Window sizes of 3 x 3 and 5x5 are selected for smoothing using sigma filtering. 
Results of smoothing obtained by 3 x 3 sigma filter are better than those of moving 
average filter. 

4.2.3 Median filtering 

Median filtering is a nonlinear signal processing technique which is generally 
used for image smoothing and has received much of the attention due to its edge 
preserving nature. Tukey (1971) has suggested it as a tool in time series analysis. 
Afterwards, it was applied in digital image processing. Median filtering can be 
one-dimensional or two-dimensional. It is performed on an image space by moving 
a window over the pixel points and replace the central pixel of the kernel with 
the median of the original grey levels within the window. The output image thus 
received is smoother than the original one. 
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The classical image smoothing or noise suppression is considered to adopt a 
low-pass filter and in many cases this is found to be a more appropriate procedure. 
However, in the case of remotely sensed satellite imagery, it could not yield the 
satisfactory results as it fails to preserve the most vulnerable feature of the image 
i.e. edges or lineaments. In certain situations median filtering is better and two of 
its main advantages are : 

(1) It preserves sharp edges while linear low-pass filtering blurs such edges. 

(2) It is very efficient for smoothing spiky noise. This property is illustrated in 
Fig. 4.1. 

One-dimensional median filtering 

The median filtering may be applied in one-dimension case by selecting a one- 
dimensional window and moving it on to the image space. Considering the size 
of the window as m with its members as xj, . . . , x^- These members are arranged 
in order and if m is odd then grey level corresponding to the middle number is 
designated as median and if m is even then median is computed as the mean of the 
two middle numbers. 

Median of a sequence is expressed as Median(xi, . . . , x^). 

A median filter of size m on a sequence Xi,i G Z and m being an odd number, 
is expressed as 

yi — Median Xi = Median (xi_„, . . . , x,-, . . . , x,+.„), i G Z. (4.5) 

where t; = (m — l)/2 and Z corresponds to the set of all natural numbers. It 
is obvious from above definition that median filter preserves edges, whereas the 
corresponding moving average filter transforms an edge into a ramp of width m as 
shown by the following equation. 

Zi — (xi_„ -j- . . . -h Xj -T . . . + Xi+„)/m, i GZ (4.6) 
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Two-dimensional median filtering 

In most of the image processing operations it is found that digital images are 
represented on a square grid by grey levels. A two-dimensional image is expressed 
as /(x,y). 

Where (x, y) lies on or subset of 

The image is considered two-dimensional with specific numbers of rows and 
columns as in the case of satellite imagery. For a filter window (A) of size u x v, 
median of an image /(x, y) is defined as 


ff(x, y) = Median{/(x, y)} = Median[/(a; -|- u, y -h v); (u, v) € A], (x, y) G Z^ 

(4.7) 

The selection of the shape and size of the filter window depends on the antic- 
ipated feature classes present in the image. Windows may have a shape like line 
segments, squares, circles, crosses, square frames and circle rings. Some of these 
windows are shown in Fig. 4.2. 

Edge preservation and noise reduction, the two important characteristics of the 
median filtering are explained below : 

(a) Edge preservation by median filter 

An edge or boundary in image processing gives an impression of an image in 
which all pixels on one side of the line have a common grey level value p and all 
points on the other side of the line have a another common value q, and p ^ q. 
Using this definition of edge, the edge preservation property of the median filter is 
presented by the following result : 

If the median filter window set A is symmetric around the origin and includes 
the origin, i.e. if 

(r, s) e A (-r, -s) G A, (4.8) 
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Fig. 4.2 Configuration of filter windows for median filtering; 

(a) Horizontal line segement (size 3); (b) Vertical line segement (size 3) 
(c) Horizontal line segement (size 5); (d) Vertical line segement (size 5) 
(e) Plus shape window (size 3); (f) Plus shape window (size 3); 

(g) Plus shape window (size 3); (h) Plus shape window (size 5). 



(0, 0) € A, (4.9) 

then the median filter preserves any edge picture (Huang, 1981). It should be noted 
here that these conditions imply that the niunber of points n in A is odd. 

Since a symmetrical window is considered for filtering and median is being cal- 
culated within that window, it is obvious that boundary of the two regions will be 
preserved. 

(b) Noise reduction by median filter 

Basic aim of median filter is to suppress the noise. Pratt (1991) has discussed in 
rather qualitative way its effect on impulse noise. Fig. 4.3 illustrates some examples 
of the operation of a median filter and a moving average filter for a discrete step 
function, ramp function, pulse function and a triangle function with a window of 
size five. It is seen from these examples that median filter has the desirable property 
of not affecting step functions or ramp functions. Pulse functions with periods less 
than one-half the window width, are suppressed. But, the peak of the triangle is 
flattened. 

Operation of the median filter can be analyzed qualitatively to a certain extent. 
The median of a product of a constant K and a sequence f{j) is expressed (Pratt, 
1991) as 

Median [K{f(j)}] = K Median {f(j)} (4.10) 

Median [K + f(j)] = K -f- Median {f(j)} (4-11) 

This result can also be generalized for two dimensional case /(cc, y) as follows: 


Median [K{f(x,y)}] = K Median 


(4.12) 

Median [K "l|j^(x, y)] = K + Median 


(4,13) 
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Fig. 4. 3 Moving average and Medlon filtering on various test signals 



Thus it is evident that salient details do not change when median is performed 
on the raw image or even on translated or transformed ( by a factor K) image. 
4.2.4 Pseudomedian filtering 

In comparison to moving average filtering, the median filtering is found to be 
more efficient as far as edge preservation and the noise suppression axe concerned. 
But its computational requirements are higher as number of operations, for every 
slide of the window, required are arranging members within the window in order, 
locating the middle number and replacing the central pixel by it. Number of com- 
putations increase if window size is increased. Attempts have been made by earlier 
workers to evolve various methods to attain high throughput using median filtering. 
Pratt (1984) evolved a filter known as the pseudomedian filter which possesses the 
properties of the median filter. 

Consider the median of a five-element sequence x\,X 2 .,xz,Xi,xz which can be 
expressed as 

MED(xi,X 2,5C3, a:4,X5) = MAX[MIN(xi,X2)®3)5MIN(xi,X2>iC4)j 

MIN(xi, X2, ^s), MIN(xa, X3, X4), 

MIN(xi, X3, X5), MIN(xi, X4, X5), 

MIN(x2, X3, X4), MIN(x2, X3, X5), 

MIN(x2, X4, X5), MIN(x3, X4, X5)] 

or 

MED(X],X2,X3,X4,X5) = MIN[MAX(Xi,X 2,X3 ),MAX(Xi,X 2,X4), 

MAX(Xi,X 2,X5), 

MAX(xi, X3, X4), MAX(xi, X3, X5), MAX(xi, X4, X5), 

MAX(x2, X3, X4), MAX(x2, X3, X5), MAX(X2, X4, X5), 

MAX(x3,X4,X5)] (4.15) 

It is found that for a sequence of size I, the median will be the MAX ( or MIN) 

of the MIN (or MAX) of all 11 /[n! ][(/ - n)! ] subsequences 
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where, 


n = (^ + l )/2 and 

n =1(1-1)... 2. 

For a specific case of a one-dimensional array of order five, psuedomedian is 
computed as follows : 

Pmed(xi,X2,X3,X4,X5) = (l/2)MAX[MIN(xi,X2,X3),MIN(x2,X3,X4), 

MIN(x3,X4,X5)] 

+ (1/2)MIN[MAX(xi,X2,X3),MAX(x2,X3,X4), 

MAX(x3,X4,X5)] (4.16) 

Considering the same array as above, ten possible subsequences of the order 

three are used for computing the median. But in the case of pseudomedian filtering 
only the three subsequences (rci, X 2 , 0 : 3 ), (x 2 , 0 : 3 , 0 : 4 ) and (xz, X 4 , X 5 ) are used. In this 
case the MIN followed by MAX in the first half of Equation 4.16 always gives the 
median of the array or a grey level which is smaller than the median, while the MAX 
followed by the MIN gives the median or a larger grey level value. The average of 
the two grey level values are then computed with a view to cancel the discrepancy, 
if any. 

Pratt (1991) on the basis of an analysis with 120 possible arrangements of the five 
elements has indicated that in only eight of the cases the pseudomedian was found 
equal to the median. At the same time pseudomedian was never found more than 
the median by classical method in any arrangement. This indicates the effective 
elimination of the noise content, specially the spikes and stripping problem. 

Let Im denote an array of elements Zi, / 2 , • • - , - The pseudomedian of the array 

is 


Pmed(l„) = ( 1 / 2 )MAXIMIN(U) + ( 1 / 2 )MINIMAX(L) (4.17) 
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Here the MAXIMIN the MINIMAX are computed as follows by consideriug a 
subsequence of the size n 
where, n = (m + l)/2 

MAXIMIN(U) = MAX[MIN(li,...,U), 

MINQa, • • - jln+i)? 

. . . , MIN(ln,-n+l, ? 1-m)] 

and 

MINIMAX(U = MIN[MAX(li,...,l„), 

MAX(l2,...,Ui), 

...,MAX(lin-n+l>---jlm)] ('^•19) 

Since in the present study two-dimensional image or array is being considered, it 

is essential to extend the one-dimensional pseudomedian concept to two-dimensional 
satellite image case in a simple way. Array of any shape such as rectangular, square 
or plus type can be chosen. While applying the MAX or MIN operations the 
resulting image is oversmoothed hence few details disappear compared to classical 
median filtering. 

The application of a plus-shape pseudomedian generally provides a much better 
result in comparison to MAX and MIN operation or the classical median filtering 
approaches. While applying the pseudomedian filtering several configurations of the 
window shape and size were tried and it was realized that the plus shape window 
yields superior results and it is also computationally much faster. 

Several windows selected for pseudomedian filtering are illustrated in the Fig. 
4.4. 

For a plus shaped window of a general size 1, the pseudomedian is computed as 
follows : 

Pmed = (l/2)MAX[MAXIMIN(pi),MAXIMIN(qj)] 

+ (l/2)MIN[MINIMAX(p,),MINIMAX(q])] (4.20) 
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Fig. 4.4 Various window configurations utilized in performing 
pseudomedian filtering 

(a) Horizontal line segment (size 3), (b) Vertical line segment (size 3), 
(c) Solid square (3x3), (d) Horizontal line segment (size 5), 

(e) Vertical line segment (size 5), (f) Plus segment (5x5), 

(g) Plus segment (size 5). 




• „ its smoothed output for moving 

4 5 Origiue^ test regio ( ) 
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where, p; and qi are grey levels for central row and central column. 

It is found that the MAXIMIN operator in one or two-dimensional form is useful 
for removing bright spikes and impulse noise, but has little or no effect on dark 
spikes or impulse noise. On the other hand, the MINIMAX operator is very good 
at removing dark spikes or impulse noise, but not the bright ones. A solution to 
this problem was suggested as to cascade the operators (Pratt et al., 1984). The 
pseudomedian operator used here removes either of the noises very efficiently. 

Out of various configurations the best results have been obtained in the case 
of plus shape window. Plus shape windows of sizes 3x3 and 5x5 are employed. 
The output obtained by these appear similar. In present case, a 3 x 3 size has been 
selected based on computational efficiency compared to a larger window. 

Based on the above theoretical background and discussions it has been inferred 
that pseudomedian filtering will yield best smoothing result. For a visual assessment 
of these filtering techniques a small test region of the study area has been selected 
as shown in Fig. 4.5(a). The following s »«oothing filters found most statable in each 
of the categories were applied on this region: (a) Moving average filter (3 x 3) 
(b) Sigma filter (3 x 3) (c) Median filter (3 x 3) (d) Pseudomedian filter (Plus shape 
of size 3) 

Outputs received after smoothing by these filters are shown in Fig. 4.5 (b), (c), 
(d), and (e). In the original image four linear features 1 (Railway line), 2 (tributary 
to river Bijul), 3 (tributory to river Bijul) and 4 (road) have been identified for 
assessing the effect of smoothing. It is evident that 1 has completely vanished in 
(b) and (c). Linear features 2 and 3 have become blurred in (b) and (c). Linear 
feature 4 is blurred in (b) but is lost in (c). In (d) and (e) images these features are 
preserved and smoothing is carried out effectively. A comparison of median filtered 
output (d) and pseudomedian filtered output (e) revealed that these linear features 
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appeared sharply in (e). 

Therefore, based on the results of this experiment, Pseudomedian filtering has 
been found most suitable for the present case of remotely sensed images. The 
smoothing output received after applying a plus shape (size 3) pseudomedian filter 
is shown in Fig. 4.6. 



Fig. 4.6 Smoothened image using a 3 x 3 plus shape pseudomedian filtering 


4.3 Edge Detection 

In the present study the main objective is to enhance the image so that 
maximum structural details are revealed. Faults, fractures, joints and other linea- 
ments can be identified easily if edges in an image are detected. A rich repertoire 
of tools are available in classical image processing which may filter out these edges 
or linear features. For the delineation of lineaments, quite often edge detection is 
carried out using Sobel operator, Laplacian operator and Roberts operator, with- 
out clearly ju.stifying the basis of their selection for the concerned applications. 
Therefore, the main objective here is to choose the most appropriate edge detection 
algorithm for lineament detection. 

In a grey level picture an edge is defined as the boundary between two 
regions of different constant grey levels. Edges, as conceived in image processing 
are generally of the following types : (Peli and Malah, 1982; Pratt et ah, 1984) 

(1) Ideal step edge 

(2) Ideal roof edge 

(3) Ideal spike edge 

(4) Noisy edges superimposed on noisy image 

Depending upon the class of image being analyzed, a variety of edge con- 
figurations exist. An illu-stration of various kinds of edges is presented in Fig. 4.7. 

Edge detection algorithms may be broadly divided into two groups;parallel 
and sequential. 
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Fig.4-7 One- dimensional ideal and real edge signals 



Parallel detection algorithms read and process all pixels of the image space 
simultaneously and output the obtained edges in parallel. In sequential detection 
algorithms, pixels are serially read and on the basis of neighborhood pixels, the 
pixels in consideration are cissigned to either an edge or background domain (Smith 
and Davis, 1975). 

Lineaments visible on satellite imagery are very useful from resource study 
point of view. Earlier, photo-interpreters job had been to identify relevant lin- 
eaments associated for various kind of application and prepare a map for further 
study. It has been felt over the years that image processing plays a significant role 
in making this process more reliable and efl&cient. With this in view several edge 
detectors have been used time and again. A large number of software packages are 
available based on common edge detection techniques. However, they often fail to 
serve the user in diverse image application requirements. 

It has been seen that in several fields including remote sensing, workers have 
been using the filtering based on the concept of computing derivatives in an image. 
Other methods are also present such as frequency domain filtering, geometrical 
image processing, optical image processing and entropy based approaches. Filtering 
based on spatial derivative and entropy may be a viable option for remotely sensed 
imagery. In the following sections a detailed discussion and use of these algorithms 
are presented. 

4.3.1 Derivative filtering 

The edge detectors used here are based on the computation of derivatives in 
a small neighborhood of central pixel. This neighborhood can be selected depending 
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on type of image to be processed. The derivatives to be computed in these edge 
detectors are of the first and second order. These are discussed in detail in the 
following sections. 

(a) First order derivative filtering 

The operators based on this concept can be grouped into two categories. 
The first one involves computation of gradients in two orthogonal directions in an 
image. The other computes a set of directional derivatives. For a few of these 
operators, considered appropriate in present study, the details are as follows : 

(i) Orthogonal gradient generation 

In orthogonal gradient based methods, gradient in a small neighborhood 
of the image is evaluated in two orthogonal directions. Using these orthogonal 
components the magnitude and slope of the gradient are determined. In the case 
of a continuous image /(a;,y), an edge can be identified by forming the continuous 
one-dimensional gradient g{x, y) along a line normal to the edge slope. The slope 
of this is at an angle a with respect to the horizontal aods. If the magnitude of 
gradient is sufficiently large, i.e., above a prespecified threshold value, an edge is 
deemed present otherwise not. The gradient along the line normal to the edge 
slope can be computed in terms of the derivatives along orthogonal axes as in the 
following equation 


gi^^y) = 


dfjx, y) 
dx 


cosa + 


df{x, y) 
dy 


sina 


(4.21) 


Since images of concern in this chapter are discrete in nature, a scheme for 
computing the edges following the above method is presented in Fig. 4.8. 

The magnitude of combined gradient is computed as follows : 
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(4.22) 


9{x.y) 


br(®,y)]^ + bc(x,y)]' 


0.5 


where, gr{x,y) is row gradient and gc{x,y) is column gradient. 

In the present study, the above equation is used for the magnitude of com- 
bined gradient. But several workers have used the following approximate formula 
for computational efficiency. 


= lyr(a;,t/)| + \gr{x,y)\ 


(4.23) 


The slope of the combined gradient with respect to the horizontal axes can 
be computed as follows : 


e(x,y)=tan-‘(?4^} (4.24) 

Ur(x,y) J 

Several operators have been proposed to accomplish the discrete domain 
differentiation approximating the above mentioned continuous differential equation 
with a view to provide a computationally efficient and better solution than the 
earlier methods. 

One simple approach of extracting edges using derivative filtering is to 
compute rmming difference of pixels along rows and columns of the image. The row 
gradient and column gradients are defined (Jensen, 1986) as row gradient 

9 rix, y) = f{x, y) - fix, y-1) (4.25) 

column gradient ; 

9 cix, y) = fix, y) - fix -t- 1, y) (4.26) 
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An Edge 





ig. 4.8 Derivative filtering for edge detection 






The row and column gradients are computed so as to ensture that gr and 
remain positive for an edge which increases in magnitude from left-to-right and 
from bottom-to-top in the discrete satellite images. 

Of the two windows of size 2x2 and 3 x 3, It has been observed that 3x3 
provides advantage of increased smoothing over 2x2. This makes the derivative 
operations less sensitive to the noise (Gonzalez and Wintz, 1977). Roberts operator 
of the size 2x2, has been used by several workers in remote sensing but was not 
used for the reason stated above. A 3 x 3 edge gradient operator has been floated 
by Prewitt (1970) which can be described by the pixel numbering convention as 
shown below 


aO 

al 

a2 

a7 

fix. y) 

a3 

a6 

a5 

a4 


The Prewitt square root edge gradient may be defined as 

= I br(®,y)]^ + [yc(a;,y)]^| 

with row gradient component 


(4.27) 


9rix,y) 


1 


k + 2 

and column gradient component 


(a2 + kas + a4) — (oo + kar + og) 


(4.28) 


9c(x,y) 


+ 2 


(oo + kai -+- ^2) — (oe + + ^4) 


(4.29) 


where A: = 1. 


82 



In the above equation the row and column gradients are normalized to 
provide unit-gain positive weighted and unit gain negative weighted averages. 

Sobel has proposed an orthogonal edge detector which differs from Prewitt 
edge detector in a way that the values of the north, south, east and west pixels are 
doubled, i.e., k = 2 (Duda and Hart, 1973). This has been done for the purpose of 
assigning equal weights to each pixel for its contribution to the gradient magnitude. 
An algorithm proposed by Prei and Chen (1977) has given weighting to north, south, 
east and west pixels hy k = \f2 so that the gradient is same for the horizontal, 
vertical, and diagonal edges. The advantage of the Prei-Chen edge detector is that 
the edge gradient is the same at the edge center for the single pixel transition 
vertical and diagonal ramp edge situation. In the case of Sobel, the gradient for 
diagonal edge is 1.06 and for Prewitt 0.94. Thus Sobel operator is very sensitive 
to diagonal edges than the Prewitt edge detector while the Prewitt edge detector 
is more sensitive to the horizontal and vertical edges than the diagonal edges. The 
Prei-Chen operator is equally efficient in the case of horizontal and vertical edges 
as well as for the diagonal edge case. 

The computation of gradient as discussed above is linear in nature. This 
involves computation of row and column gradients, which can be represented by 
simple convolution operation. 


Pr(®,2/) = f{x,y)eK{x,y) 

(4.30) 

9 c{x, y) = /(s, y) © hc{x, y) 

(4.31) 


where. 
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hri^^^y) hci^.y) are 3x3 row and column impulse response operators 
respectively and are illustrated for present image processing case in Fig. 4.9. 

It has been realized that these operators provide noisy edges in a noisy 
image environment. This problem is alleviated by resorting to a larger window 
size. To incorporate this point, softwares are developed in a general way to adopt 
any user-specified window size. These kind of operators are also called as boxcar 
operators. In present context, a window size of 3 x 3 was selected and it has 
been observed that it yields satisfactory results. When resorted to a larger window 
dimension as in the case of boxcxir filters the computation time has increased which 
makes the process highly uneconomical and at the same time some salient edges are 
left out. 

Hence the choice has been narrowed to Frei-Chen 3x3 derivative filtering. 
Thresholding 

Large number of false edges appear when an edge detector is applied. Al- 
though the data image is smoothed to suppress scattered noise from the image space 
still some eddies and other tiny edge contours appear. To make the interpretation 
easier, insignificant edges must be eliminated. This was achieved in the present 
case using selective thresholding. Threshold calculation is based on the standard 
deviation for each of the filter window. This standard deviation is updated with 
each shifting of the window. The threshold is calculated as : 

T^p + Kxa (4.32) 
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Roberts: 


Row Gradient Column Gradient 

/o 0 -1\ /-I 0 0\ 

p 1 0 0 10 

\ooo/ \ooo/ 

Prewitt : 

Row Gradient Column Gradient 



Sobel : 

Row Gradient Column Gradient 



Prei-Chen : 

Row Gradient Column Gradient 


1 

2 + \/2 


1 0 -1 \ 
\/2 0 - v ^ 

1 0 -1 / 


1 

2 + 


-1 - y /2 - 1 \ 

0 0 0 
1 1 / 


Fig. 4.9 : Impulse response operators used in orthogonal 
derivative gradient edge detection 
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where, 

T ; threshold, 

// : mean of the window pixel grey levels, 

K : coefficient by the user and 

a : standard deviation for window pixel grey levels. 

Computation of the magnitude of the gradient is followed by threshold 
computation, which is again followed by a decision making process of assigning 
the central pixel to edge domain or to background. If magnitude of the gradient 
is greater than the threshold then the central pixel is assigned to the edge pixel 
domain, otherwise it is not. 

Based on the methodologies discussed above Prie-Chen, Prewitt and Sobel 
3x3 operators are applied on the image of a test region of the study area. The 
edge output received are shown in Figs. 4.10 (a), (b), and (c) respectively. 

It is seen that a large number of noisy edge contours appear in b and c. A 
small rectangle shown by green arrow has come out as hollow in a but is filled in 
b and c. This made it clear that FVie-Chen operator has filtered the internal edges 
more accurately than Prewitt and Sobel operator. A vertical edge (3-B and lower 
part of 3-C) has been filtered very sharply in a than in other two images. An overall 
inspection suggested that edges found by applying Prie-chen operator were thin and 
sharp while other two rendered thick edges. For lineament mapping a sharp edge is 
preferred. 

Therefore, the result obtained in this experiment suggested to adopt Frei- 
Chen operator for edge detection in present context. An edge image received after 
employing this operator is shown in Fig. 4.11. 
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Fig. 4.10 Edge detection output using first order orthogonal derivative filtering on an 
image subregion of the area (a) Frei-Chen operator (b) Prewitt operator (c) 
Sobel operator 
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(ii) Directional gradient filtering 

Geological features are generally of irregular shape. Structural features are 
commonly linear or curvilinear. Often they are oriented in few preferred directions 
within the region. Directional aspect thus is very important. If prior information 
regarding the direction is available it is easier to design or adopt a suitable filter. 
Generally this information can be obtained in the field or from available geological 
maps. If access to such information is not possible then the only choice left is to 
test the filters of ail possible orientation. 

In the previous section, the orthogonal edge gradient methods have been 
discussed. Here a different approach suggested by several researchers is adopted 
which is to compute gradients in several directions. This method is also known as 
edge template gradient generation. It is accomplished by convolution of a digital 
image with a .set of directional edge templates. The edge template gradient is 
expressed as 


g{x,y) = MAX[lgi(x,y)l,...,|gi(x,y)|,...,lgn(x,y)|] (4.33) 

where, 

is the gradient in the direction which is computed by convolving an 
image with a gradient impulse response array hi(x,y). Out of all the directional 
gradients computed for each of the setting of window in the image space the edge 
direction is computed by the direction of the maximum gradient. Prewitt (1970) 
has suggested compass gradient operator which is gain-normalized. Kirsch (1971) 
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has proposed a 3 x 3 directional gradient operator using modulo 7 operation 


9{x,y) = MAXZl 


\5pi - Zqi 


(4.35) 


where, 


Pi = ai + Oj+i + ai+2 


(4.36) 


and 


9t — Oi.f3 -f- a, -4.4 + Oj+s -f- Oj+g + a, -4,7 (^-ST) 

Here, the pixel numbering convention is same as in the Section 4.3.1. Robin- 
son (1977) has proposed 3-leveI and 5-level directional template gradients. The 
3-level Robinson operator is derived after the Prewitt operator and the 5-level 
Robinson operator is derived after the Sobel operator. In this study an applica- 
tion of Prewitt’s compass gradient, Kirsch directional gradient, Robinson 3-level 
and Robinson 5-level directional gradient have been applied with an assumption 
that there is no knowledge base regarding lineament orientation available for the 
area into coasideration. Directional gradient templates are shown in Figs. 4.12 to 
4.15. These gradient templates have been utilised to produce the edge images. It is 
found that edges received in all the four images were similar except the noisy edge 
contours received in the case of Prewitt compass gradient and Robinson 3- level 
directional gradient. The edge output received by the Kirsch directional gradient 
and the Robinson 5- level direction gradient is shown in Figs. 4.16 and 4.17. Large 
number of edges appear in these images which may be utilised for lineament map- 
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ping in case of remotely sensed images in association with available ground truth 
informations such as geological and topographical maps. 

(b) Second order derivative filtering 

In the case of second order derivative the edges appear more prominently 
than first order derivative filtering (Gonzalez and Wintz, 1977). While taking the 
second derivative at boundaries or edges, it is found that extraordinary change oc- 
curs in magnitude or in sign of the derivative. Two types of second order derivative 
methods are available 

(i) Laplacian and 

(ii) Directed second order derivative. 

(i) Laplacian operator 

This impulse response operator is based on second order derivative concept 
and is rotation invariant. When applied in the continuous domain signal 

case, the Laplacian can be expressed as follows : 

hihj) = 

where. 


== 1 

di^ df 


(4.39) 
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East (/ii) 


Northeast {h^) 



North {hz) Northwest [h^) 



West {hz) Southwest {h^) 



South {h^) Southeast (/ig) 



Fig. 4.12 Directional gradient templates for the Prewitt compass 

gradient method 
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East {hi) 


Northeast (/^2) 


1 

15 





North (/13) 


Northwest (/14) 


1 

15 



-3 

0 

5 



/-3 -3 
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West (/15) 


55 


Northeast (h^) 


/■~3 -3 
~3 0 

\-3 -3 




5 

0 

-3 


5 

5 
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South (/17) 


Southeast (hs) 


1 

15 





Fig. 4.13 Directional gradient templates for the Kirsch 
gradient method 
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East (hi) 


Northeast (ht) 



North (hz) 


Northwest (/14) 




-1 

0 

1 


0 

1 
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West ih) 


Northeast (h^) 



0 

0 

0 




1 

0 

-1 


2 

1 

0 


South (hy) 


Southeast (/is) 


1 

4 



2 

0 

-2 




Fig. 4.1s Directional gradient templates for the Robinson 
5-level gradient method 
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Fig. 4.16 Edge detected image using the Kirsch gradient 






The magnitude of Laplacian is zero if f{i,j) is constant or changing 

in a linear manner. 

If the rate of change of /(i,i) is greater than linear, fiiij) exhibits a sign 
transformation at the point of inflexion of f{i,j). This point is called as the zero 
crossing of /£,(*, i) and manifests the presence of an edge. The negative sign has 
been incorporated in the Expression 4.38 so that zero crossing of fLihj) should 
invariably possess a positive slope for an edge whose amplitude increases from left- 
to-right and from bottom-to-top in an image (Pratt, 1991). This operator can easily 
be applied to the case of discrete remote sensing imagery by computing difference 
of gradients along each axis. Mathematically, this is expressed as 




— 

f{x,y + l) - f{x,y) 


t — t 

+ 

1 

4- — .,Jt 

— 

f{x-l,y) - f{x,y) 


Since this process is essentially performed by convolution it can be repre- 
sented as 


/L(x,y) = /(a;,y) oA(a:,y) 


with, 


/ 0 0 0 \ /O -1 0\ 

/i -1 2 -1 + p 2 0 

\ 0 0 0 / VO -1 0/ 

or 


h = 



-1 

4 

-1 



(4.41) 


(4.42) 


(4.43) 
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A four ncighlxir 3 x 3 Laplacian which is also normalized to unit gain is 

expressed as 

0 -1 0 \ 

-1 4 -1 I (4.44) 

0 - 10 / 

Prewitt has sxiggested an eight-neighbor Laplacian, defined by the unit gain 
normalized impulse response operator, as 




(4.45) 


This array is not separable into a smn of second derivatives as in Equation 
4.42. A eight- neighbor LapUirian impulse response operator which is separable in 
nature can be expressed as 



(4.46) 


The gain normalized separable eight neighbour Laplacian impulse response 
operator is represented as 



(4.47) 


It has been reported earlier (Heurtas and Medioni, 1986; Pratt, 1991) for 
the response of Laplacian operator in discrete image case that the zero crossings do 
not always lie at pixel sample location. As remotely sensed imagery are occasionally 
subjected to solar luminance variability hence they have ramp edges of varying 
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slopes. Due to this zero valued Laplacian response do not occur many times which 
makes detection «if efiges liiihcult. 

A simi>le approarli to detect Laplacian zero crossing in discrete domain is 
to find maximum of all positive Laplacian responses and minimum of all negative 
Laplacian resi)onses in a 3 x 3 window. If the magnitude of the difference between 
the maxima and the minima exceeds a threshold, an edge is adjudged present. The 
Laplacian operator is sensitive to the grey level jump but it is considered over- 
sensitive to the noise in the image domain (Gonzalez and Wintz, 1977). Due to 
this nature, it is less iitilized in edge detection. It is considered useful only for 
preliminary edg<' d(*tect ion. Hence, this operator is applied on the study area data 
set, only for giving a broad idea of the lineament fabric. 

A sf'parable (fight neighbor gain normalized Laplacian operator is applied 
on the original iinag«* and tin* output edge image is illustrated in Fig. 4.18. 

(ii) Directed second order derivative generation 

In above discussion of Laplacian operator, rotationally invariant second 
order derivative has been employed to ascertain the existence of an edge. The 
direction of an <'dge may be determined during zero crossing detection process. In 
this section, an alternatj* approach to first estimate the edge direction, and then 
compute th(' (uu'-dimtmsional second order derivative along the edge direction is 
presented. The directed sr*cond order derivative of a continuous domain image 
f{ij) along a liiK‘ at an angle 0 with respect to horizontal axis is given by 




v) 2 J, ^ d^f{u,v) 

dudv 




(4.48) 


Directed second order derivative is a non-linear operator. Therefore, con- 
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volution of a sniooflunj; opfrator ovor f{u,v) before differentiation process is not 
equivalent to ronvolvilu’; the directed second order derivative of f{u,v) with the 
smoothing function. 

The main problem with this filter is to have the knowledge base of the 
slope of edge' domain pixels. Tins can be solved by selecting a suitable first order 
derivative filter and determining the slope. Selection of an appropriate filter itself is 
a difficult task. Onc{‘ the probaV)le slope of the edge pixel is determined then above 
equation which is applicable only for continuous signal case is approximated for the 
discrete image rase. 

The edge dirt'ction can be determined for each of the neighborhood case but 
the problem ari.ses in approximating the continuous function /(u, u). In the absence 
of any strong methodf>logy to implement this procedure for discrete images, it was 
decided to select only the Laplacian .second order approach out of the two available 
second order derivative hit (*ring techniques for lineament detection. 

4.3,2 Edge detection based on image entropy 

Edges are characterized by sudden variation of grey level values in am image 
space. Reunotely sense<l image.s having large number of structural features, may 
be characterized Ijy a random model involving a large family of such edge pixels 
scattered in an image. Derivative filters discns.sed above are generally employed for 
detecting the edges. They generally show larger magnitude of gradient in a zone of 
sudden variation of grey h'vc'I information. Other schemes for edge detection not 
involving the derivative filtering approaches are also available in literature {Pratt, 
1978; Rosenfield and Kak, 1982). 

In this section, an approach based on the entropy of the image is employed 
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for odg<' (if'tortinn it) rrmotrly soiisfHi imagery. Entropy of the image in a predefined 
small m'ighborh'Hxl of thf global image is computed whith yields low values in the 
areas where the grey level values are going through rapid variation rather abruptly 
(Shiozaki. 19b3). E<lges are characterized by sudden variation in the grey levels 
hence they are located in the areas where the entropy is small. The edges may thus 
be filtered out from background by operating the entropy operator throughout an 
image and (ietect t he arc’as where the entropy is small. 

The entropy Imseci approach provide.s is a non-linear spatial filter. It eval- 
uates the entropy of the image grey level values in a small neighborhood centered 
at the (r. pixel. The coidiguration of pixels in a window is as follows : 

( ai a2 m 

ag Ui 

«7 «€ 

L«’t the grey levnd.s of the (r, pixel and the n neighbouring pixels is 
00,01,02, re.spectively. 

Then the entropy E of the brightness in the region is defined as follows : 


( 4 . 49 ) 

log(n + 1 ) 


where 


P. = 


(h 


( 4 . 50 ) 


The t-echnicpies of e(ig<» detection discussed in previous sections such as the 
first order edge <ietection and second order edge detection utilize the basic concept 
of rate of change of pixel gvey level values in a local neighborhood. They are not 
dependent on th(? average grey level value of that local region. When entropy of 
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the image sub-region is determined both these factors are taken into account which 
make this algorithm more rehable. However, this aspect makes it computationally 
slow. In the present study different window sizes are considered to understand the 
effect of size of neighborhood in extracting the edge pixels. Sizes of neighborhood 
considered are 3 x 3, 5 x 5 and 7x7. It is found that 3x3 neighborhood gives better 
response in comparison to other selected window sizes. The efficiency of entropy 
based method decreases with the increase in window size. In case of 7 x 7 entropy 
filter, the number of edges were less compared to 3 x 3 and 5x5. Therefor edge 
output using 3x3 window size is selected for presented work which is shown in Fig. 
4.19. 

The entropy operator is highly susceptible to the noise. Therefore the image 
must be smoothed before edge detection to eliminate noise. 
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Fig. 4.18 Edge detected image using Laplacian eight neighbour, 



Fig. 4.19 Edge detection using non-linear entropy 3x3 operator 


103 


4.4 Edge Linking 

In previous section, edge or lineament detection techniques have been dis- 
cussed and it is found that all of them are based on detection of illumination discon- 
tinuities. They detect pixels lying at the boundary of background and the object 
or between two objects. In the case of remote sensing imagery this concept is not 
applicable, because of presence of noise, break in boundary due to topographical 
or nonuniform illumination situation and vcirious other effects that induce spurious 
intensity discontinuities. 

Due to this, edge detection algorit hms must be followed by edge linking 
and other boundary detection procedures designed to assemble edge pixels into a 
more meaningful set of object boundaries. The procedure considered here is based 
on computation of edge orientation and magnitude. 

The procedure followed for edge linking or boundary linking is dependent 
on a small neighborhood similar to edge detection process. All points having similar 
characteristics are linked, forming an almost continuous edge much better than that 
achieved by classical edge detection algorithm. The two characteristics to evaluate 
the similarity of edge pixels are : 

1. The strength or the magnitude of the gradient operator used to produce 

the edge pixel 

2. The direction of the gradient. 

The strength of the gradient as computed in earUer sections is represented 
as |^[/(a:,y)]|. 

For edge linking operation, the edge detected image received from earUer 
methods is utilized. A pixel s', y' situated in a predefined neighborhood of (s, y) 
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was analyzed for edge linking as expressed in following algorithm : 
if lsl/(x.!/)]-!(l/(a:',!,')l|<r and 
|9(x,!/) - «(x',y)| < A 
then = k 

where 

T and A are the gradient and angle thresholds respectively. 

Txy is the grey level of (x, pixel. 
k denotes the grey level of the edge pixel. 

9{x, y) is the direction of the gradient vector. 

Based on the above algorithm, a point (x', y') in the predefined neighbour- 
hood of (x, y) are linked. This operation is performed by moving the predefined 
window throughout the image. An edge detected output received using Prei-Chen 
derivative filtering was processed using edge linking operation. The output image 
is shown in Fig. 4.20. 
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Fig. 4.20 Edge detected image after edge linking (The Prei-Chen image is used for 
edge linking) 
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4.5 Results and Discussion 

A detailed assessment of image processing tools for edge detection has been 
carried out in the present chapter. It is proposed to perform image smoothing op- 
eration prior to edge detection. For smoothing the satellite images, four techniques 
have been evaluated and the result is presented in Table 4.1 


Table 4.1 Results of image smoothing techniques on a selected region of the 
image 


Name of the 

technique 

Size of 

kernel 

Smoothing 

Finer 

Uneament 

details 

Computat- 
ional time 

required 

Moving average 

5x5 

over and blurred 

mostly lost 

Very slow 

filtering 

3x3 

blurred 

some of the salient 

details are lost 

Slow 

Sigma filtering 

5x5 

over and blurred 

mostly lost 

Very slow 


3x3 

Patchy appearance 

Some details 

are lost but 

better than 

moving av. filt. 

Slow 

1 

Median filtering 


Excellent 

Preserved 

Very slow 

Pseudomedian 

filtering 

3x3 

plus shape 

Excellent 

Preserved 

Fast 


Based on the above results it is decided to opt for pseudomedian filtering for 
smoothing the satellite images in order to attain excellent smoothing, preservation of 
lineament fabric and high processing speed. The smoothed image of the study area 
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using this technique is shown in Fig. 4.6. This image has been used in subsequent 
edge detection operations. 

The delineation of lineaments is carried out by detection of edges in the satellite 
imagery. Several methodologies discussed here come from derivative filtering aj>- 
proach or entropy based approach. The results obtained are summarised in Table 
4.2. 

Based on these results, it has been realised that entropy based filtering do not 
yield encouraging response in the case of remotely sensed images of Earth surface. 
The second order Laplacian operator misses some other salient edges (ramp edges 
in nature) and hence can be utilized only for preliminary study of lineaments of 
the area. The best result has been obtained in the case of first order Prei-Chen 
derivative filtering. The edge output produced using this method may be utilised 
for lineament mapping. Results of edge detection yielded by Kirsch and Robinson 
5-level operators appear impressive. These two edge outputs may be utilized in 
association with Prei-Chen output for final mapping of the lineaments. 

Edge linking has been performed on the edge output received using the Prei- 
Chen operator (Fig. 4.11) and the edge linked image is shown in Fig. 4.20. It is 
found that some discontinuties have been filled which helps the process of lineament 
extraction and visual mapping. 

In this chapter, the potentials of various digital image processing tools available 
to extract the edges and map the lineaments are presented. Beised on these tools, 
an image processing scheme for delineation of mapping of edges has been proposed, 
which is detailed in Chapter 6. 



Table 4.2 Results of edge detection techniques on a selected region 
of the image 


Technique 

Result 

First order derivative 

filtering 

(a) Orthogonal derivative filtering 


(i) Prewitt operator 

Very good detection of horizontal and 

vertical edges, thick edges obtained. 

(ii) Sobel operator 

Very good detection of inclined edges, 

large number of small edge contours formed. 

(iii) Prei-Chen operator 

Very good detection of horizontal, vertical 

and inclined edges, edges are sharper than 

Prewitt and Sobel operators. 

(b) Directional gradient method 


(i) Prewitt compass gradient 

Large number of edges detected but a large 

number of small edges. 

(ii) Kirsch gradient operator 

Very good detection of edges, less noisy 

than Prewitt compass gradient method. 

(iii) Robinson 3-level 

Very good edge detection, smaller 

directional operator 

edges found in large numbers. 

(iv) Robinson 5-level 

Very good edge detection, smaller 

directional operator 

edges found in small nmnbers. 

(c) Second order Laplacian 

Edges are detected but some of 

operator 

the salient edges are lost. 

(d) Entropy based filtering 

Only few edges detected. 
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Chapter 5 


MORPHOLOGICAL IMAGE TRANSFORMATIONS 
FOR ENHANCEMENT OF LINEAMENTS 


5.1 Introduction 

In the previous chapter, structural features of digital satellite imagery have been 
extracted using classical image processing techniques. The other group of image 
processing techniques, being discussed in this chapter, is the morphological image 
processing. In this an image is understood to be a set of natural or real numbers 
which is constituted by a number of subsets comprising numerous different struc- 
turing elements. These subsets are having a relationship which can be described in 
the notion of set theory as proposed by Matheron (1975) and Serra (1982, 1986) 
and Giardina and Dougherty (1988). 

Morphological image processing has been utilized in many areas and has gained 
popularity in recent years. Conrad in 1972, employed principles of mathematical 
morphology for studying the permeability characteristics of petroleum rocks using 
the simulation of a set of geological faults (Serra, 1982). Photographs of dust par- 
ticles had been used to understand their size distribution and shape characteristics 
(Klein, 1973). Meyer (1977) employed mathematical morphology to build a device 
for automatically performing cytological examinations. In recent past, scientists 
in the area of remote sensing have initiated studies involving morphological image 
processing for feature identification in satellite imagery. In this chapter, the pos- 
sibility of employing the concept of morphological image processing for enhancing 
structural features using remotely sensed satellite data are explored. In remote 
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sensing this concept is new and needs a thorough testing before a formal use. The 
attempt in this chapter is to develop a sequence of operations using the concept of 
mathematical morphology to achieve the same objective as in previous chapter i.e. 
extraction of geological structural details of the area. 

Morphological operations for noise removal and edge enhancement have been de- 
veloped based on the concept of mathematical morphology. The results obtained by 
the tools of morphological image processing are compared with the results obtained 
using techniques discussed in the previous chapter. 

5.2 Mathematical Morphology 

The concept of mathematical morphology was first proposed by Matheron in 
1964, who developed the concept in a realistic manner by utilizing the morphology 
or geometric shape of the constituent elements of an object. The morphological 
image processing is generally based on the analysis of binary image in terms of 
some predetermined geometric shape known as a structuring element (Giardina 
and Dougherty, 1988). The methodology involved is also dependent on the selection 
of a suitable structuring element and investigating its behaviour in the image set. 
Development of the mathematical morphology was done utilizing the binary image 
case. Only recently, this concept heis been developed for the case of multilevel 
or greyscale images. The remotely sensed images are essentially grey scale images. 
Therefore it is now possible to explore the usefulness of the concept of mathematical 
morphology for image processing of satellite digital images. 

In mathematical morphology an image /(a:, y) is considered as a set which is a 
subset of n-dimensional euclidean space %. 

Any binary image object B is represented as follows : 

B = {z: f{z) = 1, z = (x, y), {x, y) € 71^} (5.1) 
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The binary image background is a set which is defined as ; 


S'" = {z : f{z) = 0,z= (x, y) € (5.2) 

Here the function / is known as the characteristic fimction of B. If we are using 
digitized image data set which is an element of euclidean grid then the above 
equations can be represented as 

B = {{u,v) : f{u, v) = 1, (u, v) e Z^} (5.3) 

The binary image background is a set B^ which is defined as 

H*" = {(u, v) : f(u, v) = 0, (u, v) G Z^} (5.4) 

Mathematical morphology utilizes some set theoretic notions which form its foun- 
dation for morphological image processing and are summeirized in the following 
paragraphs. 

A digital image is considered an arbitrary space B. An object B of this space is 
a subset i.e., B C B. Thus an image space is constituted by infinite number of such 
image subsets. A faunily of such subsets like B constitute a set Z which is supposed 
to be less random in nature than the image space B. Some of its properties are 
listed as follows : 

1. .2 is a complete array with a partial-ordering relation called inclusion. Any 
finite set which represent an object Bi E Z m the image space has a least 
upper boimd and a greatest lower bound both belonging to Z. 

2. The lattice Z is distributive, i.e., 

Bu{CnZ) = {BuC)n{BllZ) (5.5) 
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Each of the subset like B possesses a complement within Z. 

BUB<^ = Z and J5 n 5" = $ (5.6) 

If we have two sets B and C, their difference can be expressed as 

BjC = S n C" 

In case of two sets B and C which belong to 2, the following relations are possible 

1. 5 is included in (7 [B cC] 

2. B hits C [B C : meaning BnC ^ 

3. B misses C [B C : meaning B n C = $] 

The Boolean algebra provides the framework for such morphological treatments for 
all image objects. These basic definitions are expressed pictorially in the Figure 5.1. 
Morphological transformations 

The images are used for obtaining the information of the inherent objects and the 
efforts have been to maximize information retrieval. Generally human vision fails to 
accomplish this task efficiently due to its spectral and spatial constraints (Gonzalez 
and Wintz, 1977). Digital image processing can alleviate this limitation to a larger 
extent. Morphological transformations attempt to establish spatial relationship 
between image subsets i.e., objects. In morphological image processing the basic 
set selected to identify any larger image subset of the principal image Z is very 
important. This basic set is termed as structuring element. 

In an image space Z for each of its points i a separate structuring element a(i) 
can be associated. An image can be modified to various forms by changing the 
structuring element or the transformation types. 
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Fig. 5-1 An illuslrallon of basic morphological 

(a) B] hits C ; 82 misses C ; B 

(b) Complement of C 

(c) Difference set B and C , B/C 


operations 
is included in C 
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Few basic transformations are dilation, erosion and hit and miss transfor- 
mation. Mathematically they are shown as follows : 

Dilation of an image M 

DILATE{M} : {i : a(i) j M} (5.7) 

Erosion of an image M 

ERODE{M} : {f : a{i) C M} (5.8) 

By selecting two structuring elements ai(i) and a 2 {i), the hit-or-miss transformation 
can be expressed as the difference of M eroded by ai(i) and M dilated by a^ii). 

HMT{M} : (ERODE{M}by a{) / (DILATE{M}by af) (5.9) 

Several other morphological transformations ^ are possible on these lines which may 
use union, complement, dilation, erosion, hit or miss, and various other available 
set theoretic notions. 

An image consists of wealth of information which complicates the process of 
feature identification. In morphological image processing the attempts are directed 
to manage information loss through successive transformations. These transforma- 
tions must adhere to some basic properties to yield desired quantitative and quali- 
tative results. These properties have been evolved by Serra (1986) and have become 
popular as the basic principles of mathematical morphology. These principles are 
discussed as follows : 

1. Increasing : 

If inclusion is preserved the morphological transformation is called increas- 
ing, i.e., when, 

BcC-^<S{B)c^iC) yB,C^Z (5.10) 
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2. Anti-extensivity ; 


A morphological operation is called anti-extensive when it shrinks B, i.e., 

^BeZ (5.11) 


3. Idempotence : 

After each subsequent application if the output image is unchanged by mor- 
phological transformation then such transformations are called idempotent. 

^ [if(B)] = ^{B) 'iBeZ (5.12) 

4. Homotopy : A transformation is called homotopic if it preserves the homo- 
topic sequence of an image or its subset {B). This concept of homotopy is 
illustrated in Fig. 5.2. 

The image set can be considered as a topological space. Each bounded set 
M is associated with some homotopy tree special to its kind and whose trunk 
corresponds to its background Mq (which is constituted by infinite number 
of M'^ component). The first branch of this homotopic tree is corresponding 
to Ml of M adjacent to Mq and the second branch M 2 of M next to Mi as 
illustrated in the figure. A transformation is termed £is homotopic only if it 
preserves this homotopic tree. 

5.3 Greyscale Morphological Image Processing 

Morphological image processing is a type of processing in which the spatial form 
or structure of the objects within an image are modified. Dilation and erosion are 
two fundamental morphological operations (Serra, 1982). With dilation, an object 
grows uniformly in spatial extent, while with erosion, an object shrinks uniformly. 
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Fig. 5.2 Illustration of the concept of homotopy in 
mathematical image algebra 
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The principles of mathematical morphology described in preceding sections are 
applicable for binary images and can also be extended to the greyscale image case. 
Greyscale images are described as a grey level function f{x,y) on the points of 
Euclidean 2-space. In euclidean 3-space the grey level function is considered a set of 
points [x, y, f(x, y)] imagined as a thin, imdulating and not necessarily a continuoiis 
sheet. A greyscale image is described in terms of mathematical morphology by an 
ximbra U{f) in Euclidean 3-space , where a point p with its coordinates {x,y,z) 
belongs to the xmabra if and only if z < f{x, y). In image processing the umbrae are 
very significant as they remain umbrae even after morphological operations such as 
union, intersection, dilation or erosion (Sternberg, 1986). 

Sternberg (1986) has extended the fundamental morphological operations of 
dilation, erosion, opening and closing from the binary image case to the greyscale 
images and the same are explained in following subsections. 

5.3.1 Dilation 

The dilation of a set X by structuring element Y can be expressed as : 

A:©r=UAj, (5.13) 

y&Y 

The dilation of U[a] by B is the union of translations of U[a] by points of B. At 
a given location (x, y) the grey levels d(x, y) of the dilation are determined by the 
maxima of the translated umbra U[a\. Although, U[a] is translated by all the points 
of J5, only those points on the surface of B can create maxima. The dilated greyscale 
image is computed then, as the maximum of the sum of grey levels of h(x,y), the 
function describing the surface of B, with each of the points of a(x, y), 

d(x, y) = Tnaxij[a{x -i,y-j)+ b(i,j)] (5-14) 
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5.3.2 Erosion 


The erosion of a set X by structuring element Y can be expressed as 

xey=n^-v (5-15) 

The erosion of U[a\ by JB is likewise determined as a difference of the grey levels 
of h{x, y) and points of a{x, y), the grey levels of the eroded image e(x, y) deter min ed 
by the minimum over all the points of b{xy)^ 


e(x, y) = minij[a{x -i,y- j) - b{-i, -j)] (5.16) 

The expressions for greyscale dilation and erosion bear a marked similarity to the 
convolution integral encountered frequently in digital image processing, with sums 
and differences replacing multiplication and minimum and maximum replacing sum- 
mation. Image processing through iterative morphological transformation is a pro- 
cess of selective information removal where irrelevant image content is irrecoverably 
destroyed, enhancing the contrast of the essential image features. 

5.3.3 Opening 

The opening of X by Y" is defined in terms of an erosion followed by a dilation. 
Mathematically it can be represented as: 

XK = [xeF]ey. (5.17) 

The opening of X by Y is the imion of translation of Y completely contained in X. 

5.3.4 Closing 

The closing of X by Y is a dilation followed by an erosion. Mathematically it 
can be represented as: 

X^ =[X®Y]eY. (5.18) 
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The closing of X by Y" is the complement of the translations of Y completely 
contained in 

Greyscale images opened and closed by 3-dimensional structuring elements are 
transformed in predictable and useful ways. Rotationally symmetric structuring 
elements are often employed for openin g transformations because the results do not 
depend on the orientation of the image. 

5.4 Morphological Image Smoothing 

As discussed earlier satellite images consist of numerous cover types which sup- 
press structural details of the area. This problem is solved by filtering the noise 
content of the image. All scattered tiny features less than a particular size are con- 
sidered noise as they add to the misinterpretation. Morphological methods seem to 
be most appropriate to define noise and extract it more easily without disturbing 
the basic structure of the image. In following sections two methods used for image 
smoothing are discussed. 

5.4.1 Morphological image smoothing using erosion 

Pixel set constituting noise can be characterized and identified based on their 
specific spectral non-conformity in a greyscale image ceise as well as by their size. A 
structural element of small size either binary or multilevel can be selected which can 
be used for erosion of greyscale image. In present study a 3 x 3 structural element 
as shown below is selected to erode the original image space. 

/I 1 1\ 

111 
U 1 iJ 

Here the central pixel is replaced by the minimum of all the pixel considered within 
the structuring element window convolving with the image. This way any scattered 
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noise present in the image can be eliminated or eroded. This operation is based on 
the greyscale erosion concept developed by Sternberg (1986). 

This operation provided good result in the case of binary images (Salim, 1991) 
but could not produce similar result in the case of multilevel satellite imagery. The 
output received by this method appeared blurred and hazy. Therefore, this method 
has not been used in present greyscale image case. 

5.4.2 Morphological image smoothing using neighbourhood criterion 
A morphological technique based on neighbotnrhood relation has been utilized. 
The structural element selected for defining neighbourhood relation is of plus shape. 
Identification of noisy pixels and smoothing has been carried out as shown below 
IF 


THEN 


\{axy — Oi)! < Threshold; Vi = 1 to 8 


ELSE 


dxy — 


a^y = Mean of all neighbourhood pixels except the central pixel 
The pixel numbering configuration is as follows: 


/ai 

0,2 


as 

Oxy 


\tt7 


as/ 
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For each of the window cases this operation is performed and the window is 
moved over the image space to smooth it. With this method smoothing of the 
image could be achieved without compromising the loss of minor linear features. 
The method is very fast as simple comparison and difference operations are involved. 
The degree of smoothness may be controlled by asigning suitable threshold. 

Based on this algorithm the IBS band 3 image has been processed and the 
smoothened output is shown in Fig. 5.3. 

5.5 Morphological Edge Detection 

There are several transforms available to detect edges or boundaries in an image 
based on the principles of mathematical morphology. Most of them are based on the 
operations of dilation and erosion. These transforms are generally non-reversible in 
nature. They use a sequence of operations involving dilation, erosion and difference 
for detection of edges. Since extensive computations involving complex mathe- 
matics like differentiation,' integration, inverse and division are not involved these 
techniques consume lesser computer time. These techniques use simple comparison, 
addition and subtraction processes therefore the shape of object is not disturbed. 

Bsised on the morphological transformations and selection of suitable structuring 
elements, several edge detectors are being discussed in this section. 

5.5.1 Edge detection using different structural elements 

In the case of edge or lineament detection in satellite imagery using morpholog- 
ical filtering the following methodologies involving dilation and erosion have been 
utilized. 
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Fig. 5.3 Morphological i 
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1. In this process the digital satellite image f{x, y) is first dilated using a struc- 
turing element h{i,j). The original image is subtracted from this image. The 
resultant transformed image g{x, y) reveals a large number of structural fea- 
tures or the edges. It is noticed here that edges mainly of outer boimdary of 
the objects are detected. 

Step One : 


Step Two : 


DILATE [fix, y)] = fix, y) © /t(i, j) 


(5.19) 


g(x,y)= DILATE [/(x. j,)] - /(x, ») (5.20) 

Therefore, this process is relevant for extracting the geological boundaries of 
various formations and other structural details marked by outer peripheral 
boundary. 

2. In the second method the digital satellite image fix,y) is eroded using a 
structuring element hii,j). This image is subtracted from the origiaal image 
eind the resultant image g' revealed a large number of edges. 

Step One : 


Step Two : 


ERODE [/(x,y)] = fix, y) 0 


(5.21) 


g'ix, y) = fix, y) - ERODE [fix, y)] (5.22) 
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Here, edges mainly of interior boundary of the objects are delineated emphati- 
cally and the edges at the outer boundsiry are not revealed as efficiently as in first 
methodology. Since irmer edges are present in abimdance, they add to the confusion 
while interpreting the lineaments. 

The structuring element h(i,j), indicated in the Equation 5.21, is the basic 
element influencing the flnal output by its size and shape. Several such elements 
can be designed to obtain output closer to the desired one. Some of the structuring 
element used for delineating the structural features in digital image are presented 
in Fig. 5.4. 

Edge detection based on above two methods have been performed using these 
structuring elements and the best extraction of edges were obtained in the case of 
plus shaped structuring element of size 3 (Figs. 5.5 and 5.6). 

5.5.2 Edge detection using top hat transformation 

Edge detection may also be achieved using transformation involving operations 
of closing and difference. This technique is known as top hat transformation. 

Top hat morphological transformation extracts topographical peaks or valleys 
from a digital images. It can also be used in the case of a greyscale satellite image. 

Consider an image transformation involving open and difference operations per- 
formed on a raw image /(x, y). 

g{x, y) = /(x, y) - /pp(x, y) (5.23) 

Here, fpvi'X’^y) is a transformed image of the original image f{x,y) which has 
undergone open operation using a structuring element of convex disk D. 
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Fig. 5.4 Structuring elements used in edge detection 
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Fig. 5.5 Edge image using morphological edge detection (using DILATION 

operation) 



Fig. 5.6 Edge image using morphological edge detection (using EROSION 


operation) 
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Serra (1982) has discussed the geometrical and mathematical properties of open 
and close operations in detail. They are graphically indicated in Fig. 5.7. As seen 
in the figme, the changes in shape show that open, operation using a structuring 
element mentioned above trims down the capes, peaks and saddles. The close 
operation on the other hand fills the valleys and emphasizes the isthmuses. The 
difference operation performed between the original image and the image received 
after opening is shown in Fig. 5.7(c). This image reveals the safient features such as 
the top portion of the capes and peaks after filtering. In Fig. 5.7(d), the difference 
operation performed on the closed image and the original image is shown. Here, 
it may be indicated that sudden drop in elevations, discontinuities and valleys axe 
filtered out efficiently. In the present study, for the extraction of structural features 
such as faults, folds, fractures and Uneaments, a transformation involving closing 
and difference with original image is expected to yield better result than the one 
involving opening and difference operations. 

Meyer (1977) proposed the top hat transformation in his pioneering work on 
extraction of chromatin in cell nuclei. He conducted a large number of experi- 
ments which prove the power of this algorithm in extracting topographical peaks, 
isthmuses and capes. 
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Fig. 5.7 (a) OPEN operation on one-dimensional signal 

(b) CLOSE operation on one-dimensional signal 

(c) Difference (f - fg); 

g 

(d) Difference (f - f). 
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The top hat transformation can be expressed mathematically as follows : 

Top hat (/) =/*-/! (5,24) 

where, is an image obtained after performing close operation on image /. 

The top hat transformation should be performed in all possible directions since 
the orientation of linear features in a new image is not known. 

This transformation is performed in present study using a sequence of closing 
by three structuring elements having the directions (0°, 60°, 120°) and also the com- 
plementary structuring elements in the directions (30°, 90°, 150°) followed by the 
differencing operation with the original image. The structuring elements xised in 
this study have hexagonal shape as shown in Fig. 5.8. 

The smoothened image (Fig. 5.3) has been used for detection of the edges. The 
result of transformation using structuring element in 0 ° direction on the input image 
appeared very noisy. The other directions in which the transformation was applied 
were 60°, 120°j 30°, 90° and 150°. The respective edge images for these are indicated 
in Figs. 5.9 to 5.13. 

5.5.3 Edge detection using top hat transformation and image superim- 
position 

In the preceding section, edge detection based on top hat transformation has 
been discussed. This technique has brought out edge details efficiently which is 
useful for mapping the lineaments. Lineaments such as faults, fractures, dissected 
zones and joints are delineated on the satellite imagery using certain visual keys like 
isolated linear vegetation trend, drainage pattern and moisture in comparatively 
dry region. These cover types are associated with lower grey levels. If these cover 
types are superimposed on to the edge output generated in previous subsection then 
interpretation of lineaments would be easy. 
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This object has been achieved by involving the operation of top hat transfor- 
mation and image superimposition. The top hat transformation is performed in 
the same way as shown by Equation 5.24. The process along with the modification 
suggested here is as follows: 

(1) Perform the close operation on the image; 

(2) Subtract the original image from the output image of step 1, it is same as 
Top hat transformation; 


Tophat (/) = (5.25) 

(3) Compare each pixel of top hat (/) and original smoothened output and 
replace the pixel of top hat image if there exists a minima in the original smoothened 
image. 

EDGE(/) = MIN[{Tophat(/)}, /] (5.26) 

The operations involved in this method are indicated in Fig. 5.14. 

This method using several structuring elements as shown in Fig. 5.15 has been 
attempted and the best results have been obtained in the case of 45° , 135° and plus 
shape structuring elements of size three. In the case of other structming elements 
of larger size, fewer edges were obtained and other configuration of structuring 
elements rendered noisy edge environment unsuitable for edge mapping. 

Edge images obtained using plus shaped , 45° and 135° structuring element of 
size 3 are shown in the Figs. 5.16 to 5.18 respectively. 
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Figure 5.8 Hexagonal structuring elements chosen for directional top hat 

transformation 




Fig. 5.10 Top hat transformation at 120° using hexagonal structuring element 
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Fig. 5.12 Top hat transformation at 90*^ 
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hexagonal structuring element 






Fig. 5.13 Top hat transformation at 150° using hexagonal structuring element 
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image superimposition 
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Figure 5.15 : Structural elements used for Top hat transformation and Image 

superimposition 
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Fig. 5.16 Top hat transformation and image superimposition using 3 x 

structuring element 


3 plus 
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Fig. 5.17 Top hat transformation and image superimposition using 3 x 3 45° 

structuring element 



Fig. 5.18 Top hat transformation and image superimposition using 3x3 135° 

structuring element 


139 


5.6 Results and Discussion 

Two morphological image smoothing operations have been attempted and the 
one involving neighbourhood criterion has yielded better smoothing. This method 
is faster compared to median or pseudomedian filtering. In addition to this, the 
lineament fabric is preserved and also the appearance of overall image is sharper 
than moving average and sigma filter. The results obtained are comparable to the 
median or pseudomedian filter. 

Three transformations dedicated for edge detection have been employed. The 
first one involved image differencing and dilation or erosion. Exterior boundaries 
have been extracted using dilation and image differencing whereas difference of 
original image and eroded image offers interior edges. The former transformation 
involving dilation is more useful in case of remotely sensed images. Several struc- 
turing elements have been used but the acceptable results are found with plus shape 
element (3x3). 

The second approach adopted for edge detection is Top hat transformation. 
Hexagonal structuring elements of size two are employed to enhance the edges in 
various directions like 0°, 60°, 120°, 30°, 90° and 150°. The edge output showed 
dark zones offering the impressions of relief. The discontinuities and lineaments axe 
better perceived in these edge images than edge output generated in the Chapter 
4 using digital filtering method. Edge output for various directions offer a promise 
for better delineation of edges if used jointly. 

The third approach of edge detection included image superimposition as a mod- 
ification to top hat transformation. Three structuring elements (45°, 135° and plus 
shape (3 x 3) have been utilized. If size of structuring element is increased the 
blurring effect appeared in the edge output. 

The output obtained using this process offered most encouraging edge image on 
which faults and lineaments could be easily identified. This method provided effi- 
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dent tool to distinguish between physical objects (roads, railway lines and builtup 
features) and geological structmral featxires. 
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Chapter 6 


LINEAMENT MAPPING AND ANALYSIS 


6.1 Introduction 

The image processing tools attempted in Chapters 4 and 5 have been utilized 
to evolve two lineament detection schemes. They are Digitedge Scheme and Mor- 
phedge Scheme. Based on these schemes, lineament maps are prepared and sta- 
tistical analysis of the lineaments is performed. The results are compared with 
available published work. Based on the procesed satellite data, geological map is 
prepared. 

The methodology adopted is presented. The results axe compared with existing 
maps and with Limited field verifications. A land use / land cover map using FCC 
of principal components is prepared. The ecological impact of surface coal mining 
and massive industrialization in the study area is discussed. 

6.2 Lineament Mapping 

In earlier Chapters 4 and 5 an attempt has been made to digitally enhance the 
satellite images to reveal the maximum structural details of the study area. Based 
on the methodologies developed and discussions made two schemes of Uneament 
enhancement are proposed. The first is Digitedge Scheme which is based on the 
techniques of Chapter 4. The other scheme is designated as Morphedge Scheme 
and is based on the morphological concepts described in Chapter 5. 

In Digitedge Scheme the 512 x 512 image of IRS-LISS 11 Band 3 of the study 
area is smoothed using pseudomedian filtering. This image is further utilized for 
lineament enhancements. Prei-Chen operator is applied to this image to yield the 
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edge enhanced image which is further thresholded to eliminate insig nifi cant edges. 
In order to fill the minute gaps, the operation of edge linking is performed. The 
output image obtained is presented in Fig. 4.20. This image has been utiHzed 
as the base image to map the lineaments. Robinson 5-level and 3x3 entropy 
operators are employed to generate the output images (Figs. 4.17 and 4.19). This is 
termed an auxiliary image set deemed necessary to assist the process of lineament 
identification. The geological maps and toposheets of the study area have been 
used for verification during the process of mapping. A fiow chart illustrating the 
processess involved in this scheme is shown in Fig. 6.1 and the lineament map 
prepared based on this scheme is presented in Fig. 6.2. 

The Morphedge scheme of lineament mapping is based on the methods discussed 
in Chapter 5. Here, the same image is smoothed using morphological neighborhood 
relationship. This image (Fig. 5.3) is utilized for further edge enhancement pro- 
cedures. The top hat transformations using various shape and size of structural 
elements have been performed. These transformations are categorized in three 
groups. The top hat one has used 3x3 size diagonal structural elements at 45° 
and 135°. The top hat two has utilized a 3 x 3 plus shape structural element. 
Images of this group are further modified by performing image superimposition as 
discussed in Chapter 5. This process yielded three final images which act as input 
set for mapping the lineaments. 

The top hat three transformation has used hexagonal structuring element of 
size two. As discussed eairlier this transformation is performed at 60°, 120°, 30°, 90° 
and 150°. This process yielded five images which acted as an auxiliary image set 
to assist the process of delineation of lineaments. Using the input image set and 
auxiliary image set, lineament mapping has been carried out. 
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Fig. 6.1 Digitedge Scheme for lineament enhancement and 
mapping 
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Lineament map of study area prepared using Digitedge Scheme 




As in the earlier scheme, geological maps and toposheets of the study area 
are utilized for verification of interpretations. An outline of operations involved is 
presented as a flow chart in Fig. 6.3 and lineament map prepaired using Morphedge 
scheme is illustrated in Fig. 6.4. 

6.3 Lineament Analysis 

Lineament map prepared by Digitedge scheme shows 186 lineaments as com- 
pared to 221 lineaments in the map prepared using Morphedge scheme. In these 
maps lineaments appear to be oriented in nature. Twelve zones identified indicating 
preferred orientation of lineaments are shown in the map (Fig. 6.5). These zones 
appeau: approximately similar to the geological divisions and coal mine blocks of the 
study area (Fig. 1.2). 

To confirm the preferred orientation of lineaments in different zones, statistical 
anaJysis of the lineament data hais been performed. The length of lineaments along 
with the angle from north is meatsinred from both the maps. Considering lineaments 
ais unit vectors and taking these angles, the sines and cosines of individual vectors 
have been computed and presented in a tabular form in Appendix-A. 

Since each zone comprises a set of such unit vectors the dominant directions are 
computed from the vector resultants. The X and Y coordinates of the unit vector 
having the angle 9i is computed as 

Xi = cos 9i 
Yi = sin9i. 

The vector resultants in X and Y directions are 


Xr = J2 cos 9i 

i=l 
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Morphological Smoothing 
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.6-4 Lineament map of study area prepared using Morphedge Scheme 







K' = ^sin0i. 

i=l 

With the help of these resultants the mean direction is calculated as 

a = tan — . 

Xr 

The magnitude of resultant is computed as 

R=^/xf^ 

The magnitude of resultant is standardized by dividing R by number of obser- 
vations (n) thus obtaining mean resultant length (R). 

R = R/n 

The circular variance (Si) is computed as 

52 = 1 - 

With this information the circular standard deviation ((7c) is estimated as 

Based on the mean resultant length (R) concentration pareimeter («) for each zone 
is computed from chart (Batschelet, 1965). 

The mean direction, orientation, mean resultant length, circular variance, circu- 
lar standard deviation, concentration parameter and confidence angle as computed 
for both the lineament maps are presented in Tables 6.1 and 6.2. 

(a) Test for dispersion 

The mean resultant length (R) generally ranges from 0 to 1. This produces 
a measure of dispersion. A large value of R indicates that lineaments are tightly 
bimched with small dispersion. 
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Table 6.1 : Statistical parameters for lineament analysis (Morphedge Scheme) 


ZONE 

Mean 

Direction 

(deg.) 

Mean 

resultant 

length 

Circular 

variance 

Std. 

deviation 

Cone. 

parameter 

(k) 

Confidence 
angle 
(deg) _ 

1 

-73.65 

0.95788 

0.04 

0.04 

12.76610 

2.27 

2 

89.30 

0.78291 

0.22 

0.20 

2.64613 

9.38 

3 

81.45 

0.91879 

0.08 

0.08 

6.53940 

8.26 

4 

54.60 

0.59695 

0.40 

0.39 

1.51574 

15.55 

5 

-74.09 

0.85469 

0.15 

0.14 

3.91072 

7.19 

6 

74.46 

0.62386 

0.38 

0.36 

1.60044 

13.91 

7 

40.96 

0.87015 

0.13 

0.12 

4.17703 

17.35 

8 

-83.09 

0.76877 

0.23 

0.22 

2.54686 

11.36 

9 

83.00 

0.74124 

0.26 

0.25 

2.28930 

10.67 

10 

-52.30 

0.91276 

0.09 

0.08 

5.85220 

5.84 

11 

54.32 

0.76120 

0.24 

0.23 

2.45490 

13.25 

12 

44.61 

0.65529 

0.34 

0.33 

1.78953 

9.50 


Table 6.2: Statistical parameters for lineament analysis (Digitedge Scheme) 


ZONE 

Mean 

Direction 

(deg.) 

Mean 

resultant 

length 

Circular 

variance 

Std. 

deviation 

Cone. 

parameter 

(k) _ 

Confidence 

angle 

(deg) 

1 

-73.33 

0.96137 

0.04 

0.04 

12.76610 

2.27 

2 

-69.34 

0.88185 

0.12 

0.11 

4.48876 

10.88 

3 

-87.29 

0.99797 

0.00 

0.00 

50.24210 

3.06 

4 

49.62 

0.77259 

0.23 

0.21 

2.54686 

9.63 

5 

-86.07 

0.76747 

0.23 

0.22 

2.45490 

9.33 

6 

-74.22 

0.90537 

0.09 

0.09 

5.30470 

8.72 

7 

-55.48 

0.98441 

0.02 

0.01 

25.25220 

5.75 

8 

-83.88 

0.82597 

0.17 

0.16 

3.14262 

11.25 

9 

73.47 

0.89212 

0.11 

0.10 

4.85871 

7.36 

10 

89.54 

0.68895 

0.31 

0.30 

1.95357 

12.35 

11 

65.24 

0.61215 

0.39 

0.38 

1.55738 

22.18 

12 

-55.70 

0.74622 

0.25 

0.24 

2.36930 

9.19 
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It is obvious that values are high. This shows that the lineaments in all the zones 
are tightly bunched and have very small dispersion. Contrary to mean resultant 
length the circular variance has opposite meaning. A lower circular variance signifies 
lower dispersion and tightly bunched pattern. In the present case, it is evident that 
circular variance for all the zones in both the lineament maps are very nominal 
ranging from 0.002 to 0.403 (Tables 6.1 and 6.2). This implies that the lineaments 
are having noinor dispersion. 

The lineament density (A) has been calculated from the following formula 



where 

XI 1 = total length of lineaments, and 
A = area of the zone. 

The results are shown graphically in Fig 6.6. 

It is seen that lineament density is highest in zone 3 comprising of major bound- 
ary fault. In zones 4, 8, 9, 10, 11 and 12, the lineament density is quite high as 
compared to zones 1, 2 and 5. It is interesting to note that zones 4, 8, 9, 10, 11 and 
12 are having huge coal deposits whereas later zones are barren of coal. Thus the 
coal bearing zones consist larger density of lineaments than the adjacent regions, 
(b) Test of randomness of lineaments 

Von Mises distribution is a circular equivalent to Gaussian distribution. This 
has two parameters, a mean direction (a) and a concentration parameter (/c). It is 
considered unlmodal and symmetrical about the mean direction (Davis, 1986). 
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Lineoment density 


Lineament density in various zones 



1 2 3 4 5 6 7 8 9 10 11 12 

Zones 


Fig. 6.6 Lineament density in various zones using Moiphedge and Digitedge Scheme 
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As the concentration parameter increases, the probability of having directions 
of lineament very close to mean direction increases. This specifies a trend to the 
lineament pattern. On the other hand if concentration parameter is zero then all 
directions are possible and lineaments belong to a random phenomenon. Concen- 
tration parameters (k) are listed in Tables 6.1 and 6.2. 

This distribution is utilized to test the randomness of lineaments. 

The null hypothesis is that 

Ho : K = 0 

Hi : K> 0, 

ie if concentration parameter is zero then observations may have any direction 
and have no fixed trend, on the other hand k > 0 specifies a trend or prientation. 
Prom Tables 6.1 and 6.2, it is evident that for all the zones in both the schemes, 
concentration parameter is not only greater than zero but also has high values thus 
indicating a particular trend. 

Critical values of R at 5% significance level for Rayleigh’s test as computed firom 
the chart (Mardia, 1972) are shown in Table 6.3. 

The computed mean resultant length of lineaments are found more than the 
critical value of R. Thus according to the Rayleigh’s test the null hypothesis is 
rejected and a firm conclusion of having a preferred trend of lineament is established 
for all the zones. 

(c) Test for a specified trend 

In order to test that a population of vectors is characterized by a specified mean 
direction, Stephens (1969) proposed the use of extensive charts for determining 
critical values of direction. A simple approach is to determine confidence angle 
around the mean direction of samples and find if this interval encompasses the 
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Table 6.3 Critical orientation interval encompassing lineament directions 


Zone 


Critical orientation interval 


Morphedge Scheme 

Digitedge Scheme 


Lower limit 

Upper limit 

Lower limit 

Upper limit 

1 

-78.10 

-69.20 

-68.89 

-77.78 

2 

70.91 

107.69 

-48.01 

-90.68 

3 

65.25 

97.65 

-81.29 

-93.28 

4 

24.12 

85.08 

68.49 

30.75 

5 

-88.18 

-59.99 

-67.78 

-104.37 

6 

47.20 

101.71 

-57.14 

-91.30 

7 

6.95 

74.97 

-44.22 

-66.74 

8 

-105.35 

-60.84 

-61.84 

-105.92 

9 

62.09 

103.91 

87.88 

59.05 

10 

-63.75 

-40.85 

113.74 

65.34 

11 

28.34 

80.3 

108.71 

21.77 

12 

25.98 

63.24 

-37.69 

-73.71 
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hypothetical mean direction of samples (Davis, 1986). 

Confidence angle (<Se), a parameter to evaluate the standard error of the mean 
direction is expressed in radians as 

iS'e = Ijy/n R K 

Using this parameter a probabilistic limit of population mean direction C 2 m be 
estimated. Here an assumption is made for the estimation, that errors are normally 
distributed and interval 

a ± Zff Se 

should comprise of the true population mean direction (3% oi the times. To allow 
a 95% probability to have a population mean the interval is computed as 

a±1.96 5e 


Here a and 5e are in degrees. 

The values of confidence angles and the critical intervals encompassing lineament 
directions are illustrated in Table 6.3. 

The strike direction of bedding planes in zone 1 is essentially WNW as depicted 
in geological maps. This is in conformity with the values shown in Table 6.3. 

Zone 2, comprising of granite and gneisses, is having orientations in ESE and NE 
to ENE. This is also within the computed interval of orientation from the lineament 
map using Morphedge schemes. In Digitedge scheme the trend did not match as 
the number of observations (lineaments) were very little on the map produced using 
this scheme as compared to Morphedge scheme . 

Zone 3 is mainly comprising of the major northern boundary fault between 
Barakars and metamorphics. It is having major trend in EW direction and deflects 
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ESE in the latter half (Raja Rao, 1983). This is in conformity with the inference 
made in Table 6.3. 

2x)ne 4 essentially consists of Jingurdah area. This area is marked by number 
of faults in N-S, NNE and ENE directions. Directions obtained are 30.75° - 68.49° 
for Digitedge and 24.12° - 85° for Morphedge schemes (Table 6.4). 

Zone 5 lies to the west of Jingurdah area. It has lineeunents oriented in NW 
(59.99° - 88-17°) according to Morphedge scheme and NW (67.77° - 104.36°) ac- 
cording to Digitedge scheme. 

Zones 6 and 7 have trend in NW and SE according to both the schemes of 
lineament mapping. 

Prom earlier works in field (Raja Rao, 1983) it is foimd that beds in western 
part of study zone have a general NW - SE trend. 

Zone 8 contains approximately Nigahi and Jayant blocks. The beds have NW 
and SE trend (Raja Rao, 1983). Values in Table 6.4 also indicate their mean 
directions as NW and SE. 

Zones 9 and 10 are within the Dudhidiua and Khadia blocks, and are having 
ENE trend along with critical values of direction being NE and SE (Table 6.4). 
This is supported by the field observations that general strike of beds is in NE - 
SW direction. 

Zone 11 comprises a part of Khadia block having mean direction as 65.25° (Digit- 
edge scheme) and 54.32° (Morphedge scheme). Earlier workers have found out that 
strike in eastern part of Khadia area is trending in NE direction which supports the 
above result. 

2k>ne 12 comprises Bina, Kakri and Marrak coalfields. Here the mean direction 
of lineament is found as 44.61° in Morphedge scheme. The critical Hmit of angles 
ranges from 25.98° to 63.23°. The lineament patterns are governed by the strike of 
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beds and fractures which are having a NE orientation. The Bina block is having 
NS bed strike direction which changes to N and SW in Marrak region. Within the 
Kakri block, this takes a swing to NE (Raja Rao, 1983). 

6.4 Geological Map from Digitally Processed Satellite Data 

In Chapter 3, structural mapping of the study area has been carried out using 
FCCs and LANDSAT-TM Band 7 imagery. It has been realized that for better 
delineation of structural features, a dedicated approach of image processing is re- 
quired. Based on the mathematical morphology, a scheme of enhancing structural 
features is proposed here. 

Images generated by top hat transformations in Chapter 5 have been used in 
following scheme to map out the faults, joints, fractures etc. The procedure is 
outlined with the help of a flow chart in Fig 6.7. 

The top hat images at 60°, 120°, 30°, 90° and 150° are the input images for 
delineation of faults. An auxiliary image set was utilized to confirm the fault 
and discriminate it from other linear cultural features. This set comprised the 
top hat superimposed images using plus shape structuring element and diagonal 
structuring elements at 45° and 135°. A map prepared based on these images is 
presented in Fig 6.8. As indicated in the map, a large number of faults within the 
coal mine zones could be delineated along with the structirral lineaments in slate 
and phyllite zone in northern part of the study area. Since mining is in progress, 
extensive geological mapping is conducted in field for each of the mining blocks by 
CMPDI. The maps prepared by this organization show faults, joints etc. in detailed 
manner. The details indicated in Fig 6.8 are compared with available map (D.G./ 
3535, CMPDI). It is found that most of the interpretations are correct. 
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Morphological Smoothing 



Fig. 6.7 Image processing scheme to map structural features 





























6.5 Accuracy Estimation of the Geostructural Map 


A statistical estimation of the accuracy o! the final geostructural map shown in the Fig. 
6.1S has been carried out in order to evaluate the usefulness of the proposed digital image 
processing techniques lor lineament enhancement and delineation. The final map is 
compared and analy.scd with respect to the base map of the area (Fig. 1 .2) taken as the 
reference. 

The 12 faults present on the base map have been assigned identification numbers as F| , F2 , 
Fi , ... , Fi 2 . These are also marked on the final geostructural map. The length of the 
faults and their orientation with respect to north in the base map as well as in the final map 
are listed in the Table 6.4. The percentage linear and angular accuracies of the observed 
faults on {\\c final geostructural map as estimated with reference to the base map of the 
area using the following expressions are presented in Table 6.4. 

Absolute Linear Error, = I (Lb - Lf) / Lb I 
Percentage Linear Accuracy = (1 - 8 l) * 100 
Absolute Angular Error, £0 = I (Ob - Of) / 0b I 
Percentage Angular Accuracy = (1 - £9) * 100 


where 

Lb = the length of fault in the base map, 

Lf = the length of fault in the final geostructural map, 

01, = the orientation of the fault from the north in the base map, anti 

0 ,; = the orientation of the fault from the north in the final geostructural map. 



1 able 6.4 Estimation ol accuracy of the final geostructural map 


Fault 

Bast 

Li, (m) 

Map 

Oh ( deg) 

Final Map 

Abs. 

Linear 

Error 

Abs. 

Angular 

Error 

Linear 

Aceuracy 

Angular 

Accuracy 

U (m) 

0,. (deg) 


Ee 

% 

% 

F'l 

3240 

71 

2850 

62 

0.12 

0.13 

88 

87 

F. 

2268 

102 

2250 

100.5 

0.10 

0.10 

90 

90 

F< 

6372 

92 

4950 

92 

0.22 

0.00 

78 

100 

Fa 

3780 

102 

3150 

95 

0.17 

0.07 

83 

93 

F. 

2160 

94 

2550 

91 

0.18 

0.03 

82 

97 

F/ 

2160 

2 



- 

- 

- 

- 

F/ 

1836 

173 

- 

- 

- 

- 

- 

- 

F« 

1188 

66 

1350 

66 

0.14 

0.00 

86 

100 

F, 

1088 

40 

1050 

51 

0.03 

0.28 

97 

72 

F,o* 

1188 

52 

- 

- 

- 

- 

- 

- 

F„ 

6480 

91 

7950 

100 

0.23 1 

0.10 

77 

90 


not present in the final geostructural map 


It is observed that out of the 12 faults indicated in base map, only 9 could be demarcated 
in the final map. In tlic case of fault FpJI'ig. 6.8), its orientation N 50" as delineated in the 
final map tallies with its orientation indicated in the base map. The angular accuracy thus 
works out as 100 percent. However, the length of this fault as inferred in the present work 
and indicated on the final map is 2250 m while in the base map it was shown as 1080 m in 
length. It may be noted that this fault Fp is located in Jingurdah coal mine area which has 
undergone intensive opencast mining since 1982 (when the reference base map was 
prepared) resulting in the additional exposed length of the fault. In general, the faults (Fj, 
Fs , F |2 ) located in the Jingurdah coal mine area have been infemed to have enhanced 
lengths on the final map compared to plots on the base map. 


Using the proposed morphedge scheme for lineament detection, 36 additional new faults 
have been delineated whieh were not present on the relerenee map. Field checks in 
selected areas have confirmed their presence. 




Chapter 7 

SUMMARY AND CONCLUSIONS 


The study area forms a part of the Singrauli coal field. Typical of the Gondwana 
( oal basins, the biiigrauli basin exliibits distinct tectonic setting. The northern 
boundary is a faultcMi one bringing the Precambrian rnetarnorphic formations in 
juxtaposition with the Barakar and Raniganj formations (Permian) of the lower 
Gondwana sequence. Besides tliis boundary fault, several other intrabasinal normal 
faults of varying magnitude and with a discernible pattern are recognized. These 
an'., in general, aligned alng NW - SE, NE - SW and E - W trends. According 
to Verma and Singh (1980), there is a genetic relation between the faults and the 
basement structure in the Gondwana basins. In the Singrauli coal basin, the tilt 
of i,h(^ .sculinnuits towards the boundary fault in the north is very clearly observed 
from the nature of the dips of the formations. Thus, in the study area, structural 
S(d.ting is comiikix resulting in diverse directional features. 

While choosing the digital image processing techniques for evaluation in the 
prc^sent study as tools for structural mapping, the structural setting of the coal field 
was kept in mind especially in regard to the directional trends of the linear features 
whicji n(H:essitat(^ dire.ctional and nonlinear filtering. 

In the present study, data from four bands of IRS IB - LISS- and band 7 of 
LANDSAT-TM have been utilized. Since the data belong to the two different 
satellites, preprocessing has been performed to establish compatibility. Geomet- 
ric rectification of the data has been conducted tising cpiadratic rnapiiing function. 
Image resampling of LANDSAT-TM data was carried out using nearest neighbour- 
hood technique. In order to reduce the dimensionality of dataset and to preserve 
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( h<' maxiiimiii iiilonuaUou content, principal component analysis is performed. The 
data distributed in five bands could be compressed into three components (1, 2 and 
4) containing 94.16 percent information. 

Iwo false colour composites using IRS DATA (Bands 1,2,4 and 1,3,4) are gen- 
eiatrui. In a<ldition, jjrincijxil components 1,2 and 4 were used for the generation of 
FCC. These 1 CC’s have enabled the preparation of lithological- structural map and 
th(' geornorphological map of the area. The formation boundaries, as also the EW 
trending major boundary fault along with several intraformational faults, could be 
located visually. 

S(w('ral geomori)hic units are identified such as the structural ridges in the Pre- 
cambrian terrain, mesa and butte featmes within the Barakars along with pediment 
'/OIK'S and pedijilains. Structural control of the drainage is confirmed. 

To bring out the structural featmes clearly, an extensive study of tools of digital 
image proci'ssing was undcu’taken based on wliich a Digitedge scheme for lineament 
(h'tection is projiosed. This scheme involves image smoothing, edge detection and 
edge linking operation. Pseudomedian filtering for image smoothing was found to be 
efficient in the removal of sporadic noises and in preserving the edges. Moving aver- 
age; and sigma filtering did not yield satisfactory results. Median filtering rendered 
equally good result as of pseudomedian filtering but involved more computational 
<;im(;. Edg(;s are; e;xtrae:te;el using Fre;i-Che;n ejperator. Tliis operator has detected 
edges equally well in horizontal, vertical and inclined directions. Prewitt operator 
could not detect several edges in inclined directions while Sobel operator failed to 
detect some edges in horizontal / vertical directions. Edges produced by both these 
operators were thick compared to those with Frei-Chen operator which yielded fine 
edges, facilitating mapping of the lineaments. Laplacian and entropy- based filters 
did not yield desirable result. Edge image produced by Laplacian is useful only 
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for preliminary assessment of edges. This operator fails in the case of ramp edges 
( ommonly evidenced in an imagery due to shadow effects. Entropy- based filter 
(oiild not he uschI in the {)resent study as it yielded limited edges. Edge linking 
operation is performed on edge output obtained using Frei-Chen operator for con- 
necting the minor discontinuities in edges. The final edge linked output is used for 
lineanKuit mapping. A lineament map has been prepared using this scheme. It has 
been observed that limited fault demarcation could be achieved with this scheme. 

Morphedge scheme, based on morphological transformation has been adopted. 
Morphological transformation for image smoothing and edge detection have been 
u.s(m1. Image smoothing used was based on neighbourhood criterion. Edges were 
I)reserved and noise content was reduced. The smoothed output obtained has been 
u.sed for edge detection. Top hat transformations using hexagonal structuring el- 
ements in various directions (0*’, 30°, 60°, 90°, 120° and 150°) have been used for 
extracting edges. The results obtained were found to be remarkably efficient for 
mapi)ing of lineaments and faults. 

This transformation was further modified by including image superimposition 
technique. The output obtained using this method was better suited for locating 
structural featun'.s and discriminating them with the physical objects since vegeta- 
tion, surface water bodies and moisture zones have been superimposed on top hat 
transformed image output. 

The lineaments map prepared using the Morphedge scheme has indicated higher 
d(!nsity of lineaments. In the present study, the lineament density is more in coal 
mine areas. 

A statist ical analysis of lineanuuits obtaimnl on these maps has been carried 
out. The results of test for dispersion showed that lineaments in all the zones are 
tightly bunched as evidenced from the range of circular variance (0.002 - 0.403). A 
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pr('f('n('(i (lead of liaeaaieuts has also been evidenced from the raadonmess test. 

The trend of lineaments, as obtained in the statistical analysis are in confirmity 
with the trends reported earlier (Raja Rao, 1983). It is seen that the strike direction 
of rocks in various zones is confirming within the interval of critical orientations. 
Faults delineated using the morphological transformations are indicated in a map 
{prepared during the study. 

The present work has, thus, brought out the superiority of pseudomedian fil- 
tering for image smoothing over other smoothing techniques and the suitability of 
mori)hologic:al image transformation for fault and lineament delineation compared 
to th(^ digital filtering technique. 
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APPENDIX - A 


I able A-1 : Statistical details of lineaments extracted using Morphedge Scheme 


Precambrian Formations 


ZONE # 1 (Slate and PhylHte) 

length 

(m) 

angle (0)® 

cosO 

sinO 

600 

122 

- 0.52992 

0.84805 

900 

122 

- 0.52992 

0.84805 

1650 

100 

- 0.17365 

0.98481 

1050 

100 

- 0.17365 

0.98481 

1200 

100 

- 0.17365 

0.98481 

1200 

100 

- 0.17365 

0.98481 

2400 

100 

- 0.17365 

0.98481 

1200 

100 

- 0.17365 

0.98481 

1950 

90 

0.00000 

1.00000 

900 

90 

0.00000 

1.00000 

1050 

90 

0.00000 

1.00000 

1050 

90 

0.00000 

1.00000 

1650 

110 

- 0.34202 

0.93969 

3000 

110 

- 0.34202 

0.93969 

3300 

110 

- 0.34202 

0.93969 

900 

no 

- 0.34202 

0.93969 

750 

no 

- 0.34202 

0.93969 

1200 

104 

- 0.24192 

0.97030 

1200 

104 

- 0.24192 

0.97030 

1200 

104 

- 0.24192 

0.97030 

600 

104 

- 0.24192 

0.97030 

1650 

108 

- 0.30902 

0.95106 

1200 

108 

- 0.30902 

0.95106 

1350 

108 

- 0.30902 

0.95106 

1050 

108 

- 0.30902 

0.95106 

900 

108 

- 0.30902 

0.95106 

1050 

101 

- 0.19081 

0.98163 

1650 

101 

- 0.19081 

0.98163 

2400 

101 

- 0.19081 

0.98163 

1200 

101 

- 0.19081 

0.98163 

1200 

96 

- 0.10453 

0.99452 

750 

96 

- 0.10453 

0.99452 

2550 

96 

- 0.10453 

0.99452 

1350 

96 

- 0.10453 

0.99452 

1050 

84 

0.10453 

0.99452 

1500 

84 

0.10453 

0.99452 




3000 

91 

-0.01745 

0.99985 

1200 

91 

-0.01745 

0.99985 

900 

112 

-0.37461 

0.92718 

750 

112 

-0.37461 

0.92718 

1050 

112 

-0.37461 

0.92718 

900 

120 

-0.50000 

0.86603 

2400 

94 

-0.06976 

0.99756 

450 

94 

-0.06976 

0.99756 

1350 

94 

-0.06976 

0.99756 

1500 

86 

0.06976 

0.99756 

3150 

146 

-0.82904 

0.55919 

900 

152 

-0.88295 

0.46947 

1050 

160 

-0.93969 

0.34202 

1050 

138 

-0.74314 

0.66913 

1050 

138 

-0.74314 

0.66913 

1050 

142 

-0.78801 

0.61566 


ZONE # 2 (Gneisses ) 


1650 

154 

-0.89879 

0.43837 

1500 

102 

-0.20791 

0.97815 

450 

102 

-0.20791 

0.97815 

1050 

126 

-0.58779 

0.80902 

1050 

126 

-0.58779 

0.80902 

600 

126 

-0.58779 

0.80902 

1050 

126 

-0.58779 

0.80902 

1050 

96 

-0.10453 

0.99452 

1500 

96 

-0.10453 

0.99452 

1200 

84 

0.10453 

0.99452 

900 

84 

0.10453 

0.99452 

1500 

84 

0.10453 

0.99452 

5400 

89 

0.01745 

0.99985 

600 

26 

0.89879 

0.43837 

450 

26 

0.89879 

0.43837 

450 

26 

0.89879 

0.43837 

450 

10 

0.98481 

0.17365 

900 

88 

0.03490 

0.99939 


ZONE # 3 (Boundary fault) 


2400 

— 

96 

-0.10453 

0.99452 

750 

96 

-0.10453 

0.99452 

1950 

90 

0.00000 

1.00000 

1350 

96 

-0.10453 

0.99452 

900 

96 

-0.10453 

0.99452 

1950 

90 

0.00000 

1.00000 

1050 

34 

0.82904 

0.55919 

750 

47 

0.68200 

0.73135 




Gondwana Formations 


ZONE # 4 Ranigani Formations 


900 

114 

-0.40674 

0.91355 

1050 

30 

0.86603 

0.50000 

1050 

24 

0.91355 

0.40674 

750 

24 

0.91355 

0.40674 

450 

24 

0.91355 

0.40674 

750 

32 

0.84805 

0.52992 

750 

32 

0.84805 

0.52992 

1200 

90 

0.00000 

1.00000 

750 

2 

0.99939 

0.03490 

900 

2 

0.99939 

0.03490 

750 

2 

0.99939 

0.03490 

1050 

118 

-0.46947 

0.88295 

1050 

160 

-0.93969 

0.34202 

900 

160 

-0.93969 

0.34202 

1200 

111 

-0.35837 

0.93358 


ZONE # 5 (Barren Measures and Barakars 


1650 

121 

-0.51504 

0.85717 

1050 

121 

-0.51504 

0.85717 

750 

70 

0.34202 

0.93969 

750 

58 

0.52992 

0.84805 

600 

94 

-0.06976 

0.99756 

1500 

90 

0.00000 

1.00000 

750 

156 

-0.91355 

0.40674 

600 

140 

-0.76604 

0.64279 

600 

154 

-0.89879 

0.43837 

600 

112 

-0.37461 

0.92718 

1050 

91 

-0.01745 

0.99985 

600 

92 

-0.03490 

0.99939 

750 

102 

-0.20791 

0.97815 

750 

86 

0.06976 

0.99756 

750 

110 

-0.34202 

0.93969 

1500 

32 

0.84805 

0.52992 

1200 

156 

-0.91355 

0.40674 

1050 

100 

-0.17365 

0.98481 

1500 

120 

-0.50000 

0.86603 


ZONE # 6 (Barakars) 
2250 5( 

750 11' 

750 101 

1200 1 

600 1 

1050 2^ 

900 10^ 


0.64279 

-0.40674 

-0.30902 

0.99027 

0.99027 

0.91355 

-0.20791 


0.76604 

0.91355 

0.95106 

0.13917 

0.13917 

0.40674 

0.97815 




450 

102 

-0.20791 

0.97815 

1350 

137 

-0.73135 

0.68200 

450 

137 

-0.73135 

0.68200 

1050 

137 

-0.73135 

0.68200 

1050 

106 

-0.27564 

0.96126 

900 

6 

0.99452 

0.10453 

600 

140 

-0.76604 

0.64279 

600 

10 

0.98481 

0. 17365 

900 

40 

0.76604 

0.64279 

600 

22 

0.92718 

0.37461 


ZONE #7 (Barakars) 

1050 

84 

0.10453 

0.99452 

600 

22 

0.92718 

0.37461 

300 

20 

0.93969 

0.34202 


ZONE # 8 (Barakars) 


300 

58 

0.52992 

0.84805 

450 

58 

0.52992 

0.84805 

450 

58 

0.52992 

0.84805 

600 

58 

0.52992 

0.84805 

1050 

100 

-0.17365 

0.98481 

900 

100 

-0.17365 

0.98481 

1500 

98 

-0.13917 

0.99027 

1050 

120 

-0.50000 

0.86603 

600 

142 

-0.78801 

0.61566 

450 

142 

-0.78801 

0.61566 

1500 

24 

0.91355 

0.40674 

600 

155 

-0.90631 

0.42262 

450 

140 

-0.76604 

0.64279 


ZONE # 9 (Barakars) 

750 

58 

0.52992 

0.84805 

1200 

58 

0.52992 

0.84805 

1350 

2 

0.99939 

0.03490 

600 

128 

-0.61566 

0.78801 

750 

128 

-0.61566 

0.78801 

1200 

90 

0.00000 

1.00000 

900 

91 

-0.01745 

0.99985 

1800 

90 

0.00000 

1.00000 

750 

98 

-0.13917 

0.99027 

1200 

110 

-0.34202 

0.93969 

1200 

120 

-0.50000 

0.86603 

750 

46 

0.69466 

0.71934 

750 

46 

0.69466 

0.71934 

600 

118 

-0.46947 

0.88295 

750 

42 

0.74315 

0.66913 

600 

6 

0.99452 

0.10453 



%()_ 

162 

-0.95106 

0.30902 





ZONE #10 (Barakars) 

750 

142 

-0.788()T 

0.61566 

750 

142 

-0.78801 

0.61566 

1200 

142 

-0.78801 

0.61566 

900 

152 

-0.88295 

0.46947 

1200 

152 

-0.88295 

0.46947 

600 

138 

-0.74314 

0.66913 

1050 

95 

-0.08716 

0.99620 

600 

94 

-0.06976 

0.99756 

1050 

94 

-0.06976 

0.99756 

1050 

152 

-0.88295 

0.46947 

1050 

142 

-0.78801 

0.61566 

1800 

122 

-0.52992 

0.84805 

1050 

90 

0.00000 

1.00000 

750 

90 

0.00000 

1.00000 

1200 

128 

-0.61566 

0.78801 

1200 

152 

-0.88295 

0.46947 

900 

110 

-0.34202 

0.93969 

900 

155 

-0.90631 

0.42262 


ZONE# n (Itarakars 


450 

47 

0.68200 

0.73135 

750 

47 

0.68200 

0.73135 

900 

121 

-0.51504 

0.85717 

900 

138 

-0.74314 

0.66913 

900 

38 

0.78801 

0.61566 

450 

38 

0.78801 

0.61566 

600 

38 

0.78801 

0.61566 

750 

90 

0.00000 

1.00000 

450 

10 

0.98481 

0.17365 

750 

10 

0.98481 

0.17365 



ZONE #12 (Barak ars) 

450 36 0.80902 0.58779 

300 36 0.80902 0.58779 

600 26 0.89879 0.43837 

600 26 0.89879 0.43837 

300 26 0.89879 0.43837 

300 26 0.89879 0.43837 



-0.17365 

-0.17365 

-0.88295 

-0.88295 

-0.40674 

0.93969 

0.93969 


0.98481 

0.98481 

0.46947 

0.46947 

0.91355 

0.34202 

0.34202 



600 

14 

450 

14 

450 

14 

300 

14 

1200 

166 

600 

106 

1050 

94 

1800 

160 

1050 

0 

300 

0 

450 

8 

900 

38 

600 

38 

1500 

10 

600 

10 

1650 

16 

600 

100 


0.93969 0.34202 

0.97030 0.24192 

0.97030 0.24192 

0.97030 0.24192 

0.97030 0.24192 

- 0.97030 0.24192 

- 0.27564 0.96126 

- 0.06976 0.99756 

- 0.93969 0.34202 

1.00000 0.00000 

1.00000 0.00000 

0.99027 0.13917 

0.78801 0.61566 

0.78801 0.61566 

0.98481 0.17365 

0.98481 0.17365 

0.96126 0.27564 

- 0.17365 0.98481 


Table A-2 : Statistical details of lineaments extracted using Digitedge Scheme 


Precambrian Formations 


ZONE # 1 (Slate and Phyllite) 

length 

(m) 

angle (0)® 

COS0 

sin0 

600 

109 

-0.32557 

0.94552 

900 

109 

-0.32557 

0.94552 

300 

109 

-0.32557 

0.94552 

300 

109 

-0.32557 

0.94552 

450 

109 

-0.32557 

0.94552 

450 

109 

-0.32557 

0.94552 

750 

104 

-0.24192 

0.97030 

2100 

114 

-0.40674 

0.91355 

300 

114 

-0.40674 

0.91355 

2700 

109 

-0.32557 

0.94552 

900 

109 

-0.32557 

0.94552 

600 

79 

0.19081 

0.98163 

900 

80 

0.17365 

0.98481 

1350 

104 

-0.24192 

0.97030 

1650 

104 

-0.24192 

0.97030 

450 

104 

-0.24192 

0.97030 

2400 

106 

-0.27564 

0.96126 

450 

106 

-0.27564 

0.96126 

750 

101 

-0.19081 

0.98163 

450 

101 

-0.19081 

0.98163 

900 

101 

-0.19081 

0.98163 

1050 

101 

-0.19081 

0.98163 

1200 

101 

-0.19081 

0.98163 

1200 

101 

-0.19081 

0.98163 

1500 

101 

-0.19081 

0.98163 

1200 

120 

-0.50000 

0.86603 

1050 

130 

-0.64279 

0.76604 

600 

130 

-0.64279 

0.76604 

750 

130 

-0.64279 

0.76604 

750 

130 

-0.64279 

0.76604 

750 

130 

-0.64279 

0.76604 

450 

130 

-0.64279 

0.76604 

1200 

90 

0.00000 

1.00000 

600 

90 

0.00000 

1.00000 

1350 

90 

0.00000 

1.00000 

750 

89 

0.01745 

0.99985 

750 

89 

0.01745 

0.99985 

750 

89 

0.01745 

0.99985 

1050 

89 

0.01745 

0.99985 

600 

90 

0.00000 

1.00000 

1200 

150 

-0.86603 

0.50000 




900 

118 

900 

122 

600 

122 

900 

82 

750 

156 

1500 

100 

1200 

100 

600 

99 

750 

99 

750 

99 

1200 

99 


-0.46947 

0.88295 

-0.52992 

0.84805 

-0.52992 

0.84805 

0.13917 

0.99027 

-0.91355 

0.40674 

-0.17365 

0.98481 

-0.17365 

0.98481 

-0.15643 

0.98769 

-0.15643 

0.98769 

-0.15643 

0.98769 

-0.15643 

0.98769 


ZONE # 2 (Gneisses) 



90 

0.00000 

1.00000 


91 

-0.01745 

0.99985 


98 

-0.13917 

0.99027 


142 

-0.78801 

0.61566 


144 

-0.80902 

0.58779 


140 

-0.76604 

0.64279 


70 

0.34202 

0.93969 


ZONE # 3 (Boundary fault) 


2250 

94 

-0.06976 

0.99756 

600 

89 

0.01745 

0.99985 

1350 

90 

0.00000 

1.00000 

900 

90 

0.00000 

1.00000 

750 

90 

0.00000 

1.00000 

900 

98 

-0.13917 

0.99027 

1500 

98 

-0.13917 

0.99027 


Gondwana Formations 


ZONE # 4 Raiiiganj Formations 


1350 

114 

-0.40674 

0.91355 

600 

48 

0.66913 

0.74314 

600 

50 

0.64279 

0.76604 

450 

30 

0.86603 

0.50000 

300 

30 

0.86603 

0.50000 

450 

30 

0.86603 

0.50000 

450 

30 

0.86603 

0.50000 

600 

26 

0.89879 

0.43837 

600 

26 

0.89879 

0.43837 

600 

50 

0.64279 

0.76604 

600 

40 

0.76604 

0.64279 

600 

30 

0.86603 

0.50000 

300 

30 

0.86603 

0.50000 

900 

159 

-0.93358 

0.35837 

300 

159 

-0.93358 

0.35837 

450 

20 

0.93969 

0.34202 

600 

56 

0.55919 

0.82904 

300 

86 

0.06976 

0.99756 


ZONE # 6 (Barakars) 


1200 

120 

-0.50000 

0.86603 

1350 

124 

-0.55919 

0.82904 

1200 

166 

-0.97030 

0.24192 

1200 

128 

-0.61566 

0.78801 

750 

170 

-0.98481 

0.17365 

600 

146 

-0.82904 

0.55919 

900 

12 

0.97815 

0.20791 

900 

46 

0.69466 

0.71934 

450 

54 

0.58779 

0.80902 

600 

108 

-0.30902 

0.95106 

450 

108 

-0.30902 

0.95106 

1500 

106 

-0.27564 

0.96126 

450 

106 

-0.27564 

0.96126 

450 

74 

0.27564 

0.96126 

750 

74 

0.27564 

0.96126 

1800 

69 

0.35837 

0.93358 

450 

46 

0.69466 

0.71934 

750 

46 

0.69466 

0.71934 

450 

89 

0.01745 

0.99985 

450 

90 

0.00000 

1.00000 


ZONE # 7 (Barakars) 


900 

134 

-0.69466 

0.71934 

900 

134 

-0.69466 

0.71934 




600 

450 


120 -0.50000 0.86603 
no -0.34202 0.93969 


ZONE # 8 (Barakars) 


600 

56 

0.55919 

0.82904 

600 

50 

0.64279 

0.76604 

1050 

122 

-0.52992 

0.84805 

900 

66 

0.40674 

0.91355 

750 

154 

-0.89879 

0.43837 

900 

126 

-0.58778 

0.80902 

600 

126 

-0.58778 

0.80902 

450 

98 

-0.13917 

0.99027 

450 

58 

0.52992 

0.84805 

750 

106 

-0.27564 

0.96126 


ZONE # 9 (Barakars) 


600 

86 

0.06976 

0.99756 

750 

86 

0.06976 

0.99756 

450 

58 

0.52992 

0.84805 

1050 

68 

0.37461 

0.92718 

1500 

54 

0.58779 

0.80902 

600 

59 

0.51504 

0.85717 

750 

59 

0.51504 

0.85717 

900 

59 

0.51504 

0.85717 

1200 

no 

-0.34202 

0.93969 

750 

116 

-0.43837 

0.89879 

450 

114 

-0.40674 

0.91355 

600 

14 

0.97030 

0.24192 

450 

56 

0.55919 

0.82904 

750 

88 

0.03490 

0.99939 


ZONE #10 (Barakars) 


600 

132 

-0.66913 

0.74315 

1500 

150 

-0.86603 

0.50000 

750 

150 

-0.86603 

0.50000 

750 

44 

0.71934 

0.69466 

900 

42 

0.74314 

0.66913 

600 

42 

0.74314 

0.66913 

900 

42 

. 0.74314 

0.66913 

900 

24 

0.91355 

0.40674 

750 

14 

0.97030 

0.24192 

750 

132 

-0.66913 

0.74315 

750 

86 

0.06976 

0.99756 

450 

86 

0.06976 

0.99756 

300 

86 

0.06976 

0.99756 

600 

112 

-0.37461 

0.92718 

600 

150 

-0.86603 

0.50000 

1500 

130 

-0.64279 

0.76604 



ZONE #11 (Barakars) 


1200 

144 

-0.80902 

0.58779 

750 

134 

-0.69466 

0.71934 

900 

20 

0.93969 

0.34202 

450 

10 

0.98481 

0.17365 

450 

96 

-0.10453 

0.99452 

450 

12 

0.97815 

0.20791 

450 

60 

0.50000 

0.86603 


ZONE # 12 (Barakars) 

900 

138 

-0.74314 

0.66913 

300 

138 

-0.74314 

0.66913 

300 

138 

-0.74314 

0.66913 

600 

135 

-0.70711 

0.70711 

600 

104 

-0.24192 

0.97030 

600 

104 

-0.24192 

0.97030 

600 

76 

0.24192 

0.97030 

450 

82 

0.13917 

0.99027 

600 

102 

-0.20791 

0.97815 

900 

148 

-0.84805 

0.52992 

600 

148 

-0.84805 

0.52992 

600 

166 

-0.97030 

0.24192 

750 

168 

-0.97815 

0.20791 

450 

168 

-0.97815 

0.20791 

1350 

178 

-0.99939 

0.03490 

900 

6 

0.99452 

0.10453 

1950 

174 

-0.99452 

0.10453 

1200 

118 

-0.46947 

0.88295 

900 

120 

-0.50000 

0.86603 

600 

120 

-0.50000 

0.86603 

450 

80 

0.17365 

0.98481 

600 

24 

0.91355 

0.40674 



